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ABSTRACT 

This  3-year  AFOSR-sponsored  project  “In-Situ  Adhesive  Bond  Assessment  for 
Aerospace  Structures”  (#FA9550-07-l-0016)  considered  the  optimal  application  of  ultrasonic 
guided  waves  (using  sparse  piezoelectric  transducer  arrays)  for  probing  composite  skin-to-spar 
adhesive  joints  used  in  aerospace  structural  design  for  various  forms  of  defects  such  as 
delamination  and  poorly-cured  adhesive.  An  important  part  of  the  project  was  developing 
compensation  in  assessment  algorithms  for  changing  temperature  effects  corresponding  to 
normal  aircraft  operations  (-40  to  +60  ®C).  Specifically,  this  project  considered  three 
fundamental  research  areas:  (1)  analytical  modeling  and  experimental  validation  of  wave 
propagation  that  accounts  for  temperature  effects  on  the  transducer  piezo-mechanical  properties 
over  the  flight  temperature  range  listed  above  the  transducer-host  structure  interaction,  and  the 
host  structure  wave  dispersion  properties;  (2)  further  development  of  a  Global-Local  (GL) 
method  to  simulate  the  interaction  of  the  guided  waves  with  structural  defects  in  isotropic  and 
composite  panels  that  takes  advantage  of  the  Semi- Analytical  Finite  Element  (SAFE)  approach  for 
the  “global”  simulation  that  significantly  reduces  computational  time;  and  (3)  development  of 
physics-based  (from  parts  1  and  2  above)  and  data-based  (using  a  novel  two-step  detection  and 
classification  procedure  rooted  in  outlier  analysis  and  extreme  value  statistics)  algorithms  for 
assessing  the  in-situ  state  of  a  composite  adhesive  bond  with  regards  to  delamination  and  bond  cure 
failure  modes.  All  validation  experiments  were  conducted  on  a  Carbon-Fiber  Reinforced  Plastic 
(CFRP)  [0/+45/0]s  laminate  subjected  to  the  -40  °C  to  +60  temperature  excursion. 

EXECUTIVE  SUMMARY 

This  report  summarizes  the  key  findings  in  the  three  research  areas  listed  above  and  lists  all 
publications  that  resulted  from  the  project.  This  project  funded  the  Ph.D.  degree  completion  of  two 
graduate  student  researchers  in  the  Department  of  Structural  Engineering  at  University  of  California 
San  Diego  (Ankit  Srivastava  and  Tim  Fasel),  who  were  supervised  and  mentored  by  the  PI  Michael 
Todd  and  the  Co-PI  Francesco  Lanza  di  Scalea.  A  quad  chart  summary  of  the  project  is  included  in 
Figure  1  below. 
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Figure  1.  Final  “quad  chart”  summary  of  APOSR  project  #FA9550-07- 1-00 16. 
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The  first  part  of  this  project  has  developed  a  model  to  predict  the  amplitude  spectrum  of 
pitch-catch  guided-wave  signals  in  composite  plates  probed  by  MFC  piezocomposite 
transducers.  The  model  accounts  for  the  transducer  piezo-mechanical  properties,  the  transducer- 
plate  interaction  through  shear-lag  effects,  and  the  substrate  plate  properties,  all  of  which  can  be 
parametrized  with  temperature.  Theoretical  and  experimental  results  obtained  on  a  [0/=t45/0]s 
CFRP  laminate  between  -40  and  60  °C  show  substantial  changes  in  the  predominant  so-mode 
amplitude  (as  high  as  40%  change  from  the  ambient-temperature  value).  Little  change,  instead,  is 
seen  in  the  wavelength  tuning  points  (maxima  and  zeros  in  the  response).  It  was  also  observed 
that  the  signal  amplitude  follows  two  opposite  trends  below  and  above  ambient  temperature. 

These  findings  were  applied  to  the  detection  of  bond  defects  in  a  CFRP  skin-to-spar 
adhesive  joint.  The  defects  considered  were  a  region  of  poorly-cured  adhesive  (50%  degradation 
in  stiffnesses)  and  two  different  sizes  of  disbonds.  The  study  used  an  so  ultrasonic  frequency 
(-200  kHz)  that  corresponded  to  a  wavelength  tuning  peak  for  the  MFC  transducers,  as  well  as 
to  a  mode-coupling  condition  for  the  multilayered  joint  structure.  The  strength  of  transmission 
of  this  frequency  was  used  to  monitor  bond  state  in  the  (-40  to  60  °C)  temperature  range.  It 
was  found  that  a  multivariate  Outlier  Analysis  using  8  representative  signal  features  may  be 
successful  in  detecting  the  bond  defects  within  the  temperature  excursions.  Correct  defect 
classification  was  obtained  when  the  baseline  distribution  for  the  analysis  was  partitioned  in  two 
regions,  below  ambient  temperature  (-40  °C  to  20  °C)  and  above  ambient  temperature  (20  ®C  to 
60  ®C).  In  an  on-board  system,  this  suggests  that  the  defects  examined  can  be  successfully 
detected  by  toggling  the  baseline  between  only  two  operational  conditions:  in-flight  and  storage. 
There  are  several  assumptions  in  the  present  work;  the  model  assumes  plane  waves,  thin 
transducers  and  neglects  damping  in  the  substrate  plate.  Also,  the  results  on  the  bond  defect 
detection  only  apply  to  the  specific  CFRP  skin-to-spar  joint  tested.  Different  choices  of  wave 
frequency,  Feature  Vector  dimension,  and  baseline  distribution  will  likely  be  required  for 
different  applications. 

This  second  part  of  the  project  has  demonstrated  the  use  of  SAFE  in  a  Global-Local 
framework  which  allows  extracting  information  on  guided  waves  scattered  by  defects  in  isotropic 
and  layered  anisotropic  plates  for  quantitative  structural  health  monitoring.  The  advantage  of  SAFE 
is  its  ability  to  handle  the  complex  geometry  of  layered  composites  in  a  computationally  efficient 
manner.  The  method  was  used  in  cases  of  notch-like  defects  in  aluminum  plates  and  delamination- 
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like  defects  in  layered  composites.  Both  reflected  and  transmitted  spectra  of  propagating  modes 
were  predicted.  It  should  be  noted  that  two  different  discretization  are  done  in  the  solution  of  the 
problem.  These  affect  the  accuracy  of  the  solution  in  different  limiting  regimes.  For  example,  the 
thickness  profile  of  the  out  of  plane  displacement  component  of  the  Sq  mode  in  an  Aluminum  plate 

varies  from  being  approximately  linear  at  an  fd  value  of  0,5  to  a  half  wavelength  sinusoid  at 
fd  =  3.0 ,  This  indicates  that  the  thickness  discretization  of  the  Global  part  in  Global-Local  should 
be  increased  to  adequately  model  the  increasing  nonlinear  thickness  variation  of  the  modeshape  as 
the  frequency  is  increased.  Higher  order  modes  exhibit  a  more  nonlinear  modeshape  compared  to 
lower  order  modes.  The  mode  ,  for  instance,  has  more  wavelengths  in  its  displacement  thickness 

profile  at  fd  =  5  than  has  at  the  same  fd  value.  Therefore,  in  general,  a  higher  order  mode 

needs  a  finer  thickness  discretization  at  the  same  frequency  than  a  lower  order  mode  in  order  to 
better  approximate  the  true  modeshape. 

The  local  discretization  for  the  Finite  Element  part  of  Global-Local  should  also  vary  with 
frequency  and  the  modes  under  consideration.  At  any  given  frequency,  modes  with  shorter 
wavelengths  require  a  finer  FE  discretization  for  good  approximation.  Modes  with  low  phase 
velocities  in  the  dispersion  curve,  therefore,  require  more  finite  elements  per  unit  propagation 
distance  for  an  accurate  representation  in  the  summation  solution.  For  isotropic  structures  like  an 
aluminum  plate,  since  all  the  modes  asymptotically  approach  a  constant  phase  velocity  (  Sq  and  Aq 

converge  to  the  Rayleigh  surface  wave  velocity;  higher  order  modes  converge  to  shear  wave 
velocity),  the  FE  discretization  should  be  made  finer  as  higher  frequencies  are  considered.  It  was 
also  shown  that  the  wave  time  histories  may  be  extracted  from  these  models  for  an  arbitrary 
excitation.  The  results  shown  are  limited  to  a  2-D  study.  The  extension  of  the  work  to  3-D  cases  and 
to  full,  complete  adhesive  bond  assessment  should  be  the  subject  of  future  work. 

The  third  part  of  the  project  developed  a  data-driven  approach  for  in-situ  bond  state 
characterization  to  complement  the  physics-based  modeling  and  characterization  of  the  first  two 
parts.  Of  particular  interest  was  robustness  of  the  algorithm  to  the  widely-varying  operational 
temperature  changes  noted  previously  that  aircraft  experience  during  flight  and  ground  operations. 
Attention  was  focused  on  a  combined  unsupervised  and  supervised  learning  scheme,  where  the 
difference  between  the  two  is  that  examples  of  delaminated  bonds  or  bonds  with  improper  cure 
were  obtained  to  support  the  supervised  step.  Chaotically-modulated  guided  waves  were  created 
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and  optimized  (in  a  detection  sense)  and  used  as  probes  to  perform  the  assessment  by  building  both 
time-  and  state-space  domain  models  and  using  statistical  modeling  for  performing  damage 
classification  under  Type  I/ll  error  control.  The  use  of  chaotic  modulation  is  motivated  by  the 
extreme  sensitivity  to  its  subsequent  filtering  from  scattering  of  bond  defect  changes;  rather  than 
focusing  on  specific  mode  energy  or  mode  conversion  using  single-tone  bursts,  the  rich  phase 
structure  of  chaos  was  shown  to  be  a  sensitive  probe  in  this  application.  Multiple  chaotic  ultrasonic 
excitation  formats  were  explored,  including  short-time  chaotic  wave  packets  and  long-time  chaotic 
bulk  insonification,  in  which  the  diffuse,  reverberant  wave  field  was  examined  to  identify  structural 
changes.  This  method  of  insonification,  in  addition  to  enhanced  pattern  recognition  techniques, 
allows  this  damage  detection  scheme  to  be  employed  on  complex  structural  geometries  (such  as 
joints  or  bonds)  with  which  standard  ultrasonic-based  NDE/SHM  methods  ceinnot  be  used. 

In  this  approach,  auto-regressive  (AR)  predictive  models  (linear  time  series  predictors)  are 
used  to  generate  expected  input/output  relationships  among  transducer  pairs  in  an  ultrasonic  array. 
Both  the  model  coefficients  and  the  baseline  model’s  predictive  power  (characterized  through  some 
measure  of  prediction  error)  are  used  simultaneously  to  accomplish  both  detection  and 
characterization  of  the  defects.  The  primary  feature  is  a  vector  consistency  criterion  (VCC)  that 
essentially  looks  at  the  relative  orthogonality  of  any  two  vectors  constructed  from  a  list  of  AR 
coefficients.  The  use  of  the  VCC  metric,  with  appropriate  extreme  value  statistics  models  as 
described  in  detail  below,  allowed  for  detection/classification  rates  approaching  95%  for  various 
delamination  sizes  tested  on  the  CFRP  testbed  described  above.  Robustness  to  temperature 
variability  was  accomplished  by  parameterizing  the  AR  coefficients  as  a  function  of  temperature 
and  using  minimum-distance  fitting  to  choose  the  appropriate  temperature  normalization  in-situ  so 
that  such  effects  are  wrapped  into  the  assessment  procedure.  This  approach  resulted  in 
detection/classification  rates  being  maintained  at  90%  across  the  operating  temperature  range, 

1.  INTRODUCTION 

Many  Air  Force  platforms,  e.g.,  C/KC-135,  B-52H,  C-5A,  and  C-130,  are  being  asked  to 
provide  service  past  their  originally  intended  design  lives.  Condition-based  monitoring  (CBM) 
strategies  can  allow  maintaining  such  assets  in  readiness  while  reducing  life  cycle  costs  and 
preserving  life  safety.  Adhesively-bonded  structural  joint  assemblies,  including  lap-shear,  tear 
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strap-to-skin  joints,  and  composite-to-metal  Joints  such  as  boron  epoxy  repair  patches  for 
aluminum  panels,  are  among  the  most  critical  components  of  these  and  other  aircraft  structures. 
Wear,  fatigue,  and  corrosion  (all  often  exacerbated  by  widely  fluctuating  operating 
environments),  and  even  poor  manufacturing  practice  make  such  joint  assemblies  a  common 
failure  point.  Assessing  the  structural  state  of  adhesive  bonds  is  critical  to  ensuring  the  safe 
performance  of  the  aircraft. 

Ultrasonic  testing  remains  the  most  widely  used  technique  for  bond  testing  in  the 
aerospace  industry  (Guyott  et  al.  1986).  The  popular  Fokker  bond  tester,  for  example,  is  based 
on  ultrasonic  principles.  Conventional  ultrasonic  bond  inspection  uses  normally-incident  bulk 
wave  spectroscopy  done  in  pre-flight  or  post- flight  tests  (Figure  la).  This  method  consists  of 
measuring  the  frequency-dependent  reflection  coefficients  of  ultrasonic  waves  propagating 
through  the  bonded  structure  in  the  Megahertz  frequency  range  (Guyott  and  Cawley  1986).  The 
main  drawback  of  this  technique  is  the  limited  inspection  range  that  makes  the  method  too 
“local”  to  be  practical  for  continuous  monitoring.  At  the  same  time,  “global”  structural 
monitoring  techniques,  such  as  those  based  on  vibration  data,  work  at  frequencies  often  too  low 
for  the  scale  of  typical  bond  defects,  particularly  when  the  target  flaws  are  at  early  growth  stages. 

As  an  alternative.  Guided  Ultrasonic  Waves  (GUWs)  operating  in  the  tens  to  hundreds  of 
kilohertz  range  are  being  increasingly  used  for  in-situ  monitoring  purposes.  By  exploiting  the 
waveguide  geometry  of  the  structure,  GUWs  provide  the  large  area  coverage  of  the  “global” 
methods  while  retaining  the  defect  detection  sensitivity  of  the  “local”  methods  (Figure  lb). 
GUWs  provide  a  number  of  features  that  are  naturally  sensitive  to  bond  state,  including  wave 
velocity,  amplitude,  frequency  content,  and  dynamic  state  representation  (the  latter  of  which  will 
be  a  major  component  of  this  project).  GUWs  also  offer  the  unique  possibility  of  maximizing 
the  sensitivity  to  different  types  of  bond  conditions  by  exciting  the  appropriate  cross-sectional 
mode  shapes  through  the  selection  of  the  corresponding  frequency-velocity  (dispersive) 
combinations  or  by  exploiting  novel  features  from  user-tailored  GUWs  (such  as  chaotic  GUWs). 

Work  has  been  performed  in  the  past  on  the  effects  of  temperature  on  guided-wave 
measurements  in  composite  components,  although  no  previous  study  has  examined  the  case  of 
adhesive  joints.  The  general  conclusion  is  that  temperature  variations  induce  changes  in 
amplitude  and  phase  of  guided-wave  signals.  This  was  found  as  applied  to  the  cases  of  sandwich 
panels  at  the  low  temperature  of  -90  (Blaise  and  Chang  2001),  aluminum  plates  at  the 
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moderately  high  temperatures  of  35  to  85  (Lee  et  al.  2003,  Chambers  et  al.  2006),  and 
aluminum  plates  at  extremely  high  temperatures  up  to  230  (Schulz  et  al.  2003).  Other  studies 
on  temperature  effects  have  been  conducted  in  the  context  of  reverberating  guided-wave  signals 
(acousto-ultrasonics).  These  included  ranges  of  35  to  85  (Konstantinidis  et  al.  2006),  and  5 
°C  to  40  (Lu  and  Michaels  2005,  Michaels  and  Michaels  2005).  Very  recently,  theoretical 
models  incorporating  temperature  effects  were  developed  to  predict  guided-wave  pitch-catch 
signals  recorded  by  either  monolithic  Pb(Zr-Ti)03-PZT  patches  (Raghavan  and  Cesnik  2007, 
Lanza  di  Scalea  and  Salamone  2008)  or  Macro-Fiber  Composite-MFC  patches  (Lanza  di  Scalea 
and  Salamone  2008).  The  application  of  these  theoretical  models  was  limited  to  isotropic 
aluminum  plates. 

The  present  work  examined  guided-wave  monitoring  of  composite  panels  subjected  to 
temperature  fluctuations  of  normal  aircraft  operations  (-40  °C  to  +60  °C).  The  flexible  MFC-type 
PI  patches  (Wilkie  et  al.  2000,  Sodano  et  al.  2004)  are  considered  as  wave  transducers.  A  model 
is  first  proposed  to  predict  pitch-catch  wave  signals  in  composite  plates  under  changing 
temperature.  The  results  of  the  model  are  then  used  to  guide  experiments  on  a  CFRP  skin-to-spar 
joint  with  simulated  bond  defects  (poorly-cured  adhesive  and  two  sizes  of  disbonds).  It  is  shown 
that  a  multivariate  statistical  Outlier  Analysis,  performed  on  the  amplitude  of  the  pitch-catch 
waves  propagating  through  the  joint,  is  able  to  detect  the  bond  defects  despite  the  temperature 
excursion  once  the  statistical  baseline  distribution  is  partitioned  in  as  few  as  two  portions, 
namely  below  ambient  temperature  (-40  °C  to  +20  °C)  and  above  ambient  temperature  (+20  °C  to 
+60  °C). 

Beyond  just  detection  is  the  goal  of  quantification  and  classification.  By  studying  the 
characteristics  of  the  reflected  and  transmitted  waves,  quantitative  information  on  the  defects  can  be 
obtained.  A  hybrid  formulation  was  used  in  this  work  wherein  the  finite  element  method  is 
employed  to  model  small  regions  near  the  defect  whereas  regions  away  from  the  defect  are  modeled 
using  a  suitable  set  of  wave  functions.  Goetschel  et  al.  (1982)  developed  a  global-local  finite 
element  formulation  for  modeling  axisymmetric  scattering  of  a  steady,  compressive,  incident  elastic 
wave  in  a  homogeneous,  isotropic  host  medium  with  an  axisymmetric  inclusion.  The  hybrid  method 
was  also  used  in  to  model  defects  in  a  plate  and  a  cylinder,  respectively  (Al-Nassar  et  al.  1991, 
Rattanawangcharoen  et  al.  1997).  In  the  area  of  SHM  of  aircraft  components,  the  method  has  been 
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applied  to  model  wave  interaction  with  defected  lap-shear  joints  (Chang  and  Mai  1995),  as  well  as 
notches  and  rivet-hole  cracks  in  plates  (Mai  and  Chang  1999,  2000). 

Many  aircraft  components  are  complex  in  either  their  geometry  (e.g.  multilayers,  tapered 
thickness,  etc.)  or  their  material  properties  (e.g.  anisotropic).  In  this  case  theoretical  wave  solutions 
for  the  global  portion  are  either  nonexistent  or  hard  to  determine.  The  Semi  Analytical  Finite 
Element  (SAFE)  method  can  help  handle  these  cases  because  of  its  ability  to  extract  modal 
solutions  of  complex  structures  in  a  computationally  efficient  manner  (Hayashi  et  al.  2003,  Bartoli 
et  al.  2006).  Sabra  et  al.  (2008)  demonstrated  the  application  of  the  SAEE-aided  hybrid  formulation 
to  the  detection  of  holes  in  aluminum  plates.  The  present  work  extends  the  global- local  approach  to 
model  notches  in  aluminum  plates  and  delamination-like  defects  in  composite  panels  and  suggests 
how  this  could  apply  to  delaminations. 

This  previous  active  approach  inherently  relies  upon  characteristics  of  specific  GUW 
modes  and  their  properties,  and  the  GUWs’  interaction  with  the  structure.  When  coupled  with 
SAFE  models,  very  specific  features  may  be  identified  that  can  indicate  specific  flaws. 
However,  structural  geometries  such  as  highly  curved  regions,  bolted  interfaces,  etc.,  pose 
significant  modeling  and  testing  challenges  for  GUWs,  resulting  in  difficulty  in  classifying  size 
or  type  of  damage  due  to  the  inherent  simplicity  of  the  actively  imparted  excitation  signal,  which 
is  usually  a  modulated,  narrow-band,  short-wave  pulse.  Instead,  some  researchers  have  attempted 
to  employ  bulk  insonification,  where  an  ultrasonic  source  is  excited,  and  the  resultant  long-time, 
or  diffuse,  wave  field  is  examined  to  identify  structural  changes  (Michaels  and  Michaels,  2005). 
This  method  was  investigated  as  well  to  augment  the  standard  guided  wave  method  for  structures 
with  complex  boundary  conditions  or  geometries  that  make  tracking  and  analysis  of  a  single 
propagating  mode  difficult  or  impossible.  This  approach  relies  on  data-driven  (typically,  time 
series)  models  and  change  detection  approaches  rooted  in  statistical  pattern  recognition. 

Given  these  varied  approaches,  one  might  envision  a  possible  general  taxonomy  for 
damage  detection  strategies  in  the  ultrasonic  domain  as  shown  in  Figure  2.  Global-Local, 
including  SAFE,  modeling  (physics-based)  is  critical  for  the  middle  column,  and  data-based 
modeling  is  critical  for  the  right-hand  column. 
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Figure  2.  Ultrasonic  domain  diagnostic  taxonomy. 


The  primary  challenges  facing  GUW-based  CBM  data  interrogation  methods  are  twofold:  (i) 
extracting,  in  an  appropriately  efficient  time  frame,  features  from  the  detected  GUWs  that  are 
sensitive  to  the  defects  of  interest,  and  (ii)  ensuring  that  these  features  are  reasonably  insensitive 
to  operational,  testing,  and  environmental  variability  observed  in  platform  service,  or  using 
appropriate  statistical  modeling  to  decouple  (or  at  least  minimize  correlation  of)  the  feature 
sensitivities  to  damage  and  to  this  variability.  This  project  proposed  three  research  areas 
(summarized  in  the  Executive  Summary  Section  above)  that  addressed  these  challenges.  The 
following  three  sections  will  now  describe  the  research  findings  addressing  these  challenges  in 
greater  detail,  followed  by  a  full  listing  of  publications  resulting  from  this  project,  and  concluding 
with  cited  references. 
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2.  HEALTH  MONITORING  OF  COxMPOSITE  PANELS  UNDER  CHANGING 
TEMPERATURE 

2.1  Pitch-Catch  Guided  Wave  Model  With  Shear-Lag 

The  problem  examined,  Figure  3a,  is  that  of  an  MFC  actuator  on  a  T700/5208  CFRP 
composite  plate  with  a  lay-up  of  [0/±45/0]s,  exciting  ultrasonic  guided  waves  which  are  detected 
by  a  similar  piezoelectric  patch  sensor  (pitch-catch).  The  problem  is  studied  in  2-D,  and 
assuming  plane  waves. 


(b) 


(c) 
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Figure  3.  (a)  Pitch-catch  Lamb  wave  monitoring  by  MFC  patches;  (b)  composite  MFC  type-Pl 
patch;  (c)  equivalent  model  for  MFC  type-Pl  sensor  patch. 
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In  a  previous  work  (Lanza  di  Scalea  and  Salamone  2008),  a  pitch-catch  model  was 
developed  for  PZT  and  MFC  pitch-catch  systems  applied  to  isotropic  plates.  This  formulation 
explicitly  accounted  for  the  effects  of  temperature  on  the  transducer  piezo-mechanical  constants, 
the  transducer/plate  interaction  (through  shear-lag),  and  the  plate  dispersive  properties.  The 
model  is  here  revised  and  applied  to  the  case  of  the  composite  laminate  at  hand.  This  is  done  by 
considering  an  “equivalent”  isotropic  plate  whose  so  dispersion  properties  match  those  of  the 
composite  laminate  (determined  by  a  Semi-Analytical  Finite  Element  -  SAFE  analysis)  in  the 
frequency  range  of  interest.  The  focus  is  on  the  somode  since  this  was  previously  determined  to 
be  preferentially  excited  by  the  MFC  transducers  with  a  setup  similar  to  that  of  Figure  1  (Lanza 
di  Scalea  et  al.  2007). 

2.1.1.  Shear  Stress  at  Actuator-Plate  Interface 

Shown  in  Figure  3b,  the  MFC-type  PI  transducer  is  composed  by  PZT  fibers  that  are 
unidirectionally  aligned,  embedded  into  an  epoxy  matrix,  and  sandwiched  between  two  sets  of 
inter-digitated  electrodes.  The  strain  generated  by  an  MFC  actuator  poled  along  the  in-plane 
(fiber)  direction  1  (electric  field  components  £2=  £*3=  0)  and  excited  by  a  voltage  Fappi  is: 


^  act  1 1  ^  \ 


appl 


(1) 


where  a  is  the  MFC  inter-electrode  spacing  and  d\  \  is  the  relevant  piezoelectric  coefficient. 

Assuming  the  surface-bonded  actuator  system  of  Figure  4a,  the  non-dimensional  shear- 
lag  parameter,  F,  from  Crawley  and  de  Luis  (1987)  is  given  by: 


t  bondta 


4^  +  a 
4^ 


(2) 


where  Ghond  >  tbond  ^^e  the  shear  modulus  and  thickness  of  the  adhesive  layer;  yLi  >  /«•<  >  tad 
the  in-plane  Young’s  modulus  (fiber  direction),  length  and  thickness  of  the  MFC  actuator;  the 
factor  a  is  1  or  3  for  extensional  or  bending  actuation,  respectively;  and  'F  is  the  plate-actuator 
stiffness  ratio  given  by: 
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(3) 


 E  plate  t  plate 


yLio 


It  is  well  known  that  the  shear-lag  effect  increases  with  increasing  compliance  and  thickness  of 
the  bond  layer.  Perfect  bonding  (no  shear-lag)  is  achieved  for  a  thin  and  stiff  bond  layer  where 
r  o<5 . 


(a) 


(b) 


The  spatial  distribution  of  the  in-plane  shear  stress  applied  by  the  actuator  into  the  plate 
proposed  by  Crawley  and  de  Luis  (1987),r(.r) ,  can  be  transformed  into  the  spatial  frequency 

(wavenumber)  domain,  T(k) ,  once  the  finite  length  of  the  actuator,  ,  is  taken  into  account 
(Lanza  di  Scalea  and  Salamone  2008).  The  spatial  Fourier  transform  of  the  actuated  shear  stress 
under  shear-lag  is: 


T(k)=  -2  / 


G  bond  I  a 


rw  (4  r'  +  A:'/L 


kL 


tanh(r) 


cos 


kL 


—  2sin 


kL 


(4) 


where  i  is  the  imaginary  unit  and  k  the  wavenumber.  For  the  perfect-bond  case  ( F  — >  o°  ), 
Equation  (4)  retrieves  the  “wavelength  tuning”  conditions  identified  by  Giurgiutiu  (2005),  i.e. 
zero  actuation  when  k  2  =  nn  for  «  =  0,1,2..  (i.e.  actuator  length  equals  an  integer  multiple 

of  the  wavelength  X),  and  local  maximum  actuation  for  k  2  =  K{2n-\)l  2  (i.e.  actuator 
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length  equals  an  odd  multiple  of  the  half  wavelength  /1/2),  The  zero  and  maximum  actuation 
points  are  only  slightly  shifted  from  these  ideal  conditions  for  practical  shear-lag  cases. 

The  magnitude  of  the  Fourier  transform  |t(A:)|  from  Equation  (4)  is  plotted  versus  the 

non-dimensional  wavenumber  in  Figure  5.  This  figure  considers  a  baseline  MFC  type  PI 
actuator  GPa,  =  0.2mm,  (^jj=404xl0“^^mA^),  excited  by  a  unity  voltage  ( 

f^app/=^  V),  bonded  by  an  adhesive  layer  (/^^j  =  0.01  mm)  on  an  equivalent  isotropic  plate  ( 
Ep/afe  =  30  GPa ,  =  1 .067  mm ),  under  extensional  actuation  for  so  (cx=  1).  In  the  figure  the 

“wavelength  tuning”  behavior  is  evident. 


Figure  5.  Magnitude  of  spatial  Fourier  Transform  of  in-plane  shear  stress  at  the  actuator/plate 
interface  under  shear-lag  conditions:  (a)  changing  bond  layer  stiffness,  (b)  changing  actuator 
stiffness. 


Figure  5a  examines  the  effect  of  changing  the  bond  layer  stiffness.  |t(A:)|  is  plotted  for  a 

nominal  bond  stiffness  of  Gw=  2.8  GPa  (corresponding  to  a  shear  lag  parameter  F  =  209.8), 
and  for  two  extreme  cases  of  a  very  stiff  adhesive  (Gw“28GPaor  F  =  658.7)  and  a  very 
compliant  adhesive  ( G w”  0.28  GPa  or  F  =  66.3 ).  The  plots  indicate  that  the  predominant  effect 

of  a  changing  adhesive  stiffness  is  a  slight  shift  in  zeroes  and  maxima  of  the  r(/:)  curve  (see 


inset  in  Figure  5a),  i.e.  the  wavelength  tuning  points,  with  virtually  no  effect  on  the  magnitude  of 
jr(/r)|,  i.e.  the  net  force  transferred  to  the  plate. 
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Figure  5b  examines  the  effect  of  changing  the  actuator  stiffness.  |t(^)|  is  plotted  for  a 

constant  Qi^„j=2.%G?a ,  and  with  changing  actuator  stiffness  from  a  nominal  y^,=  30GPa  ( 
E plate!  1 ).  to  a  very  stiff  actuator  (  yL=  90GPa  or  EpiatelYlct^  0.33)  and  a  very  compliant 

actuator  (  3GPa  or  EpiatelYlct~  It  is  shown  that  the  actuator  stiffness  influences  highly 
the  magnitude  of  the  stress  transfer,  with  |t(A:)|  substantially  decreasing  with  lower  actuator 
stiffness. 

Under  temperature  variations  both  Gbond  Y^t  vary.  Shown  later,  the  temperature 
range  of  interest  in  this  study  (-40  °C  to  +60  ®C)  will  have  a  substantial  influence  on  the  actuation 
magnitude,  with  little  change  in  the  wavelength  tuning  points. 

2.1.2.  Actuated  Lamb  Strain  in  Plate 

The  Lamb  wave  strain  field  generated  by  the  actuator  is  found  from  the  mode  expansion 
formulation  by  Giurgiutiu  (2005)  coupled  with  the  shear  stress  Fourier  Transform  of  Equation 
(4).  The  actuated  Lamb  strain  at  the  surface  of  the  plate  in  the  wave  propagation  direction,  jc,  is: 


£jj.(X,/‘)j  ^  /  plaif  i  6 


i{kx-(at) 


(5) 


where  the  actuated  wave  amplitude  is: 


T\ 


2g 


plate 


Nsikl 


^Dsik  )l^k 


Symm.  modes 


Antisymm.  modes 


(6) 


where  are  the  eigenvalues  (wavenumbers)  for  the  Lamb  symmetric  modes  (solutions  of 

Ds  =  0)5  are  the  eigenvalues  for  the  Lamb  anti-symmetric  modes  (solutions  of  Da  =  ^)- 

and  the  remaining  terms  are  well-known  expressions  of  Rayleigh-Lamb  theory  (e.g.  Rose  1999). 
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In  Equation  (6)  the  shear  stress  expression  is  decomposed  into  the  following  constant  term,  n , 
and  ^-dependent  term, 


2  Gr  .  lacf  E 


Ti  = 


bond  Cact 
t  bond 

,  tanh(r) 

klac,  - = - COS 


X2^k)  =  - 


-  2sin 


(7) 


T700/5208  CFRP  Composite  Lamina 

£•,=134  GPa  £2=9  GPa  G, 2=4.2  GPa  v, 2=0.26  G23=3.2  GPa  p=1530  Kg/m^  /  =  0.13  mm 

Equivalent  Isotropic  Plate 

£/?/o/e=23.7  GPa  Gp/fl/^=13  GPa 

Table  1 .  Ambient-temperature  properties  of  CFRP  lamina  and  equivalent  isotropic  plate. 

Equation  (6)  requires  the  elastic  properties  of  the  substrate  plate,  including  its  k{(0) 
dispersion  curves.  The  plate  of  interest  in  this  study  is  a  [0/±45/0]s,  T700/5208  CFRP  composite 
laminate,  with  wave  propagation  direction  in  the  90  degree  direction  of  the  top  ply  (see  Figure 
6).  The  laminate  can  be  studied  by  considering  an  “equivalent”  homogeneous  isotropic  plate 
whose  So  dispersion  properties  match  those  of  the  composite  plate  in  the  frequency  range  of 
interest.  Similar  equivalent  isotropic  models  of  composite  plates  were  considered  in  the  past  for 
guided  wave  studies  (Sohn  et  al.  2003,  Kim  et  al.  2007).  Phase  velocity  dispersion  curves,  shown 
in  Figure  4,  were  first  obtained  for  the  subject  composite  plate  using  a  Semi-Analytical  Finite 
Element  approach  (Bartoli  et  al.  2006).  The  elastic  constants  assumed  in  the  model  for 
T700/5208  are  reported  in  Table  1.  The  same  curves  were  also  obtained  for  the  “equivalent” 
isotropic  plate  using  conventional  Rayleigh-Lamb  characteristic  equations.  The  two  independent 
engineering  constants  Epiate  and  Vpiate  assumed  for  the  “equivalent”  plate  were  iteratively  changed 
until  a  satisfactory  match  was  obtained  for  the  so  mode  of  interest  (see  Figure  6).  Density  and 
thickness  of  the  “equivalent”  isotropic  and  the  composite  plates  were  kept  identical.  The  figure 
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shows  that  a  good  so  match  does  not  translate  into  a  match  of  the  antisymmetric  ao  mode,  and 
increasingly  so  with  increasing  frequency.  This  is  because  the  phase  velocity  of  the 
antisymmetric  mode,  in  the  non-dispersive  branch,  depends  upon  the  out-of-plane  shear  modulus 
that  is  typically  small  in  composite  plates  because  mainly  influenced  by  the  matrix  properties. 
Hence,  an  equivalent  isotropic  model  cannot,  in  general,  capture  the  dispersion  properties  of  both 
symmetric  and  antisymmetric  modes  of  an  anisotropic  composite  laminate. 


Freq.  [kHz] 

Figure  6.  Dispersion  curves  of  [0/±45/0]s  T700/5208  CFRP  composite  laminate  along 
propagation  direction  x  (SAFE  solution)  and  equivalent  isotropic  plate  (Rayleigh-Lamb 
solution). 
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2. 1.3.  Sensor  Response 

The  shear-lag  model  for  the  plate  -  MFC  sensor  system  is  shown  in  Figure  4b.  The  strain 
transferred  from  the  plate,  Spiaie,  into  the  sensor,  ,  can  be  found  by  using  the  study  by 

Crawley  and  de  Luis  (1987)  for  a  piezoelectric  actuator-plate  case,  and  modifying  the  boundary 
conditions  to  reflect  the  plate-sensor  case.  The  following  boundary  conditions  are  applied  at  the 
ends  of  the  sensor: 


e Plate  =  Ee  (Lamb  strain)  and  =0  at  x  =  ±  /„„  /  2  (8) 

These  lead  to  the  following  solution  for  the  sensor  strain,  ,  as  a  function  of  the  incoming 
Lamb  strain  in  the  plate,  ,  and  the  shear-lag  parameters  'F  ,  T  and  a. 


^  vp  ^ 


Esens  Ex 
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cosh 
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^  o  ^ 

y  I  sens  J 


coshr 


(9) 


The  voltage  response  of  the  PI -type  MFC  sensor  can  be  derived  by  considering  the  equivalent 
circuit  of  Figure  Ic.  Here  the  MFC  is  reduced  to  an  equivalent  PZT  with  thickness  equal  to  the 
MFC  electrode  spacing,  a,  and  length  equal  to  the  MFC  net  electrode  length,  (p  4),  where  p  is 
the  number  of  electrodes  on  a  single  face  of  the  MFC.  The  net  electrode  length  is  thus  a  fraction 
of  the  full  MFC  length,  p  le=P  Imfc,  (here  l^fc  =lsens)  where  /J  is  a  reduction  factor  depending  on 
the  design  of  the  MFC.  The  equivalent  model  also  requires  a  rotation  of  the  sensor  axes  (7,  3)  as 
shown  in  the  figure.  This  model  yields  the  following  voltage  response  of  the  Pl-type  MFC 
sensor  as  a  function  of  the  in-plane  strain  in  the  sensor,  £  : 

/  Sffx: 

2 

V  .UFC  “  S  MFC  J  Escm  ^  1 0) 
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with  sensitivity  term: 


a  {^d\\Y\^  ^/i2  Vi2  ^2 ) 

S\fhc~~^.  f  ~  7^  7  I  E  2  FTf 

^ MFC  |_^ii(l “ V12V21) “ (c/ii/i  ^^dudnVnYi  ^dnYi)] 

where  ( yf ,  Vij )  are  the  in-plane  elastic  constants  of  the  MFC  sensor  and  (  c/jj ,  c/j2 ,  )  are  the 

nonzero  piezoelectric  coefficients  and  dielectric  permittivity  term  of  the  MFC  sensor. 

2.1.4.  Complete  Pitch-Catch  Actuator-Sensor  Response 

The  pitch-catch  solution  for  a  pair  of  MFC  transducer  patches  under  shear-lag  conditions. 
Figure  3,  can  now  be  found  by  substituting  the  strain  expression  of  Equation  (9)  into  the  sensor 
response  of  Equation  (10).  After  algebra,  the  voltage  generated  by  an  MFC  actuator-sensor  pair 
is: 


where  it  was  possible  to  isolate  the  two  terms  responsible  for  the  perfect-bond  response, 
^perf_hond^  ^nd  for  the  shear-lag  response,  ,  respectively  given  by: 


R 


perf  _  bond 


k 


R 


.shear— lag 


^Isens^ 


tanh(r)cos 


(13) 

(14) 


From  Equation  (12)  the  pitch-catch  response  emerges  as  a  combination  of  the  perfect  bond 
response  and  the  shear-lag  response,  whose  effect  is  to  decrease  the  overall  signal  amplitude 
(notice  minus  sign  ). 
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2.2.  Temperature  Effects 


The  expression  in  Equation  (12)  allows  isolating  the  temperature-dependence  of  the 
properties  of  actuator,  sensor,  transducer/plate  bond  layers,  and  plate  substrate.  These  properties 
are  examined  in  the  following  sections  in  the  range  of  normal  aircraft  operations,  -40  °C  to  +60 
®C.  Unless  otherwise  specified,  a  linear  dependence  of  each  property  on  temperature  is  assumed: 

P{T)=  M  (15) 

where  P  represents  one  of  the  properties,  T  is  the  generic  temperature.  To  is  the  ambient 
temperature  (20  ),  and  dPIdT  'is  the  sensitivity  to  temperature.  The  values  considered  are 

summarized  in  Table  2  and  discussed  in  detail  below. 


MFC  transducer  properties 

E^(ro)=4.7  GPa  v,(ro)=0.45  EX7o)=70  GPa  v/'(To)=0.31 

Jar  =-0.028  GPa/°c  avjar  =0.009/ “c  a£^/a7’=o.i6GParc  av^/ar=-o.oi3 

PC 


c/ii(7’o)=404xlO-'^m/Volt 
^dy^|^T  =  \  m/(Volt  °C) 


<5^12(7'o)=-167x10‘'^ 

nvVolt 

dd.JdT  =-0.5 
m/(Volt  °C) 


eii(7’o)=15xlO-' 

Farad/m 

acii/ar  =0.043x10" 
Farad/(m  °C)  -  below 
20  °C 

aeii/ar  =0.14x10'" 
Farad/(m  °C)  —  above 
20  °C 


Transducer/plate  bond  layer  properties 
Gao„<X7’o)=3.4  GPa  aG^^/a7’=-0.013  GPa/°C 


“Equivalent”  substrate  plate  properties 
Eptoe(7o)=23.7  GPa  Gp/a,e(ro)=13  GPa 
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=-0.014  GPa/°C 


dGpi^JdT  =-0.008  GPa/°C 


Table  2.  Temperature-dependence  of  the  properties  of  MFC  transducers,  transducer/plate  bond 
layer,  and  substrate  plate. 

2.2.1.  Temperature  Effects  on  Actuator  and  Sensor  Properties 

The  temperature-dependent  properties  of  the  MFC  actuator  and  sensor  patches 
influencing  the  response  in  Equation  (12)  are:  Young’s  moduli  Y\  ^Yi  from  Equations  (2),  (3), 
and  (1 1);  Poisson’s  ratios  v,2,V2|  from  Equation  (1 1);  piezoelectric  coefficients  from 

Equations  (1)  and  (11);  dielectric  permittivity  terms  from  Equation  (11);  length  and  thickness 

dimensions.  All  of  these  terms  are  considered  in  the  following. 

Conventional  rules  of  mixture  were  used  to  calculate  the  effective  Young’s  moduli  and 
Poisson’s  ratios  of  the  MFC  patches  from  the  constituent  properties  of  fibers  (£/,  v^)  and 

matrix  (£■»,,  v^)  through  the  fiber  volume  fraction,  Vf  (Jones  1975): 


Yx=E,Vf  +  E„{\-Vf),  yI  = 


v,2  =  K/V^  +  (l-Kpv„ 


(16) 


The  temperature  dependence  of  Young’s  modulus  and  Poisson’s  ratios  for  PZT-5A  fibers  were 
modeled  according  to  Equation  (15)  using  E^{'fo)  ~  70  GPa,  v^{Jq)  =  0.31,  and  temperature 

sensitivities  of  dE^ldT=  0.16  GPa/°C  and  dv^ldT  =-0.013/”C  (Sherrit  et  al.  1999).  The 

matrix  properties  were  obtained  experimentally  from  longitudinal  and  shear  ultrasonic  time-of- 
flight  measurements  conducted  on  epoxy  coupons  in  an  environmental  chamber  with 
temperature  control.  The  measurements  were  averaged  over  two  separate  coupons  and  three 
temperature  cycles.  The  following  matrix  properties  were  obtained:  EJLTq)  -  4.7  CPa/^C, 

v„(7’o)  “  0-45  with  sensitivities  of  d  E ^1  dT  =-0.01%  GPa/°C  and  9v„/37’  =  0.009/”C.  A 
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constant  volume  fraction  f/^  =  0.4  was  assumed.  The  temperature  sensitivity  of  the  MFC 

properties  V21  then  calculated  from  Equation  (16). 

Regarding  the  piezoelectric  coefficients,  the  following  values  were  used:  (JuiTo)- 
404x10'^^  mA^olt,  ^/^(ro)  =  “167x10"*^  nVVolt,  as  obtained  from  manufacturer’s  data; 
sensitivities  dc/ii/dT=l  m/(Volt  °C)  from  Hooker  (1998)  and  dc/12 1 3T=-0.5  m/(Volt  °C)  from 
Lee  and  Saravanos  (1998)  (notice  that  dn  of  the  MFC  patch  is  equivalent  to  c/31  of  the 
monolithic  PZT  reported  in  this  reference). 

The  dielectric  permittivity  term  was  considered  bilinear  to  best  fit  the  data  by  Lee  and 
Saravanos  (1998)  {cw  of  the  MFC  patch  is  equivalent  to  ^33  of  the  monolithic  PZT  in  this 

reference).  The  ambient  value  was  Farad/m,  with  sensitivities  de^^ldT 

=0.043x10'^  Farad/(m  °C)  below  ambient  temperature,  and  3  e,j/3T  =0.14x10"^  Farad/(m  ®C) 
above  ambient  temperature. 

Patch  length  and  thickness  changes  were  modeled  assuming  a  coefficient  of  thermal 
expansion  of  5x10'^  m/(m  ®C)  along  the  fibers  and  15x10'^  m/(m  ®C)  across  the  fibers  (Williams 
et  al.  2004). 

2.2.2.  Temperature  Effects  on  Transducer/Plate  Bond  Laver  Properties 

Shear  modulus  Gbond  and  thickness  tbond  ^  both  present  in  Equations  (2)  and  (7),  influence 
the  voltage  expression  in  Equation  (12).  Ambient-cure  epoxy  adhesive  (Loctite  from  Henkel 
Corporation,  Rocky  Hill,  CT)  was  used  in  this  study.  The  temperature  dependence  of  Ghond  was 
measured  by  ultrasonic  shear-wave  velocity  measurements  on  epoxy  coupons  tested  in  the 
autoclave  between  -40  °C  and  -^60  °C.  By  averaging  the  measurements  over  two  separate 
coupons  and  three  temperature  cycles,  the  following  regression  line  was  obtained  with  a  0.96 
regression  coefficient: 

G  w(  T){GPa]  =  3 .2  -  0.0 1 3  C]  (17) 

Changes  in  thickness  of  the  bond  due  to  thermal  expansion  were  considered  negligible  and 
ignored. 
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2.2.3.  Temperature  Effects  on  Substrate  Plate 

Properties  of  the  substrate  plate  affecting  the  pitch-catch  response  include  Young’s 
modulus  Ephie  >  thickness  ,  and  shear  modulus  Gpiate  Equations  (3)  and  (6).  In  addition, 

these  properties  affect  the  dispersion  relations  k{()S). 

In  the  subject  CFRP  laminate,  the  temperature  range  of  -40  °C  to  +60  °C  affects  primarily 
the  properties  of  the  polymeric  matrix.  The  following  temperature  degradation  factor  proposed 
by  Chamis  (1984)  for  the  stiffness  of  resin  matrices  is  used; 


F  = 


^  'j'  _  7^  ^ 

T  -T 
V  ^0 


1/2 


(18) 


where  is  the  glass  transition  temperature  at  the  reference  temperature  .  The  elastic 
properties  of  the  composite  plies  were  then  calculated  as  (Chamis,  1984): 


Ex  =  EfVf  +  F_^Ej„,  E2  =  - 


EE 


Gn 


FG 


'’|2=E„V„+F/V 


(19) 


These  temperature-dependent  values  were  inserted  into  the  SAFE  model  to  predict  the  so 
dispersion  curves  for  the  subject  [0/±45/0]s  laminate.  Finally,  the  properties  Epiate  and  Gpiaie^^ 

the  “equivalent”  isotropic  plate  were  found  by  matching  the  so  dispersion  curves  at  lOdeg 
temperature  steps  from  -40  °C  to  +60  °C.  This  procedure  led  to  the  following  properties  of  the 
“equivalent”  plate:  EpiatJiTo)~^^G  GPa,  G piateiTi^~^^  GPa,  and  temperature  sensitivities  of 

dEp/aie  /  37’=  -0.014  GPaAC  and  dGphie  /  37’  =-0.008  GPa/°C.  The  density  of  the  plate  was  set  to 

Ppiaie~^^^^  kg/m^.  The  ambient-temperature  elastic  properties  for  the  composite  and  the 

equivalent  isotropic  plate  are  summarized  in  Table  1.  Figure  7  shows  the  dispersion  curves  k(cd) 
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computed  by  SAPE  for  So  and  ao  in  the  1.067  mm-thick,  [0/±45/0]s  T700/5208  CFRP  laminate  at 
the  temperatures  of  -40  °C,  +20  °C,  and  +60  °C,  and  for  waves  propagating  at  90  deg  from  the 
top  ply.  As  expected,  a  temperature  increase  causes  a  decrease  in  phase  velocity  o)lk  for  both 
modes.  It  can  also  be  seen  that,  for  the  temperature  range  of  interest,  the  shifts  in  k{oS)  values  are 
quite  small.  Hence  little  change  should  be  expected  on  the  wavelength  tuning  points  of  the  pitch- 
catch  response  spectrum. 


Figure  7.  Dispersion  curves  of  [0/±45/0]s  T700/5208  CFRP  composite  laminate  at  the 
temperatures  of  -40°C,  +20°C  and  +60°C. 

2.3.  Experiments 

The  model  of  Equation  (12)  was  compared  to  pitch-catcH  measurements  taken  on  the 
1.067  mm-thick,  [0/±45/0]s  CFRP  laminate  (Figure  8a).  The  MFC  patches  were  type  PI 
rectangular  transducers  (Smart  Materials  Corporation,  Saratoga,  FL)  with  an  active  length  of  28 
mm,  other  properties  listed  in  Table  3,  and  installed  at  a  distance  of  75  mm  from  one  another. 
The  transducers  were  bonded  to  the  plate  using  the  ambient-cure  Loctite  epoxy  adhesive.  The 
measurements  were  taken  in  an  environmental  chamber  (Figure  6b)  between  -40  °C  and  +60  °C 
at  10  °C  intervals.  Thermocouples  were  attached  to  the  plates  to  verify  that  thermal  equilibrium 
was  achieved  at  each  temperature  step.  A  3.5-cycle,  Hanning -modulated  toneburst  was  used  as 
the  actuating  signal.  The  excitation  was  swept  between  100  kHz  and  300  kHz.  The  Root-Mean 
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Square  of  the  measured  so  mode  was  computed  at  each  excitation  frequency  to  quantify  signal 
strength. 


Composite  MFC-Pl  patch 

yf  =30  GPa 

y!  =7.65  GPa  V/2=0.3 1  dn  =  404x  1 0'*" 

r//^=-168xl0'‘" 

nvV 

mA^ 

e,,=\5x\0-'’ 

0.5  a  =  0. 1  mm  lacr^28  mm 

/ac/=  0.2  mm 

Farad/m 

Transducer-plate  adhesive  layer 

Gbond-  2.8  GPa  tbond-  0.01 

mm 

Table  3.  Ambient-temperature  properties  of  composite  MFC  patch  and  of  patch-plate  adhesive 
layer. 


Figure  8.  Experimental  tests:  (a)  CFRP  laminate  with  MFC  actuator-sensor  pair;  (b) 
environmental  chamber  with  temperature  control. 


Experimental  and  theoretical  results  are  plotted  in  Figure  9.  The  curves  are  plotted 

separately  for  the  above-ambient  temperatures  ( +20 °C<T  <  +60 °C ),  and  the  below-ambient 

temperatures  (-40°C  <T<  +20  °C).  The  plots  are  normalized  to  the  relative  maxima;  since  the 
model  neglected  wave  damping  in  the  plate,  measured  and  predicted  absolute  amplitudes  cannot 
be  compared.  The  general  trend  is  consistent  between  theory  and  experiment.  The  wavelength 
tuning  points  match  well.  The  mismatch  in  response  magnitude  could  be  caused  by  the  damping 
effects  in  the  plate  and  the  plane  wave  assumption  of  the  model.  The  results  show  that  the  main 
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role  of  temperature  is  a  change  in  the  pitch-catch  response  amplitude.  The  changes  in  the  220 
kHz  response  peak  height,  relative  to  the  ambient  temperature  value,  are  on  the  order  of  30%  at 
+60  °C  and  on  the  order  of  40%  at  -40  °C. 

In  addition,  both  theory  and  experiment  show  the  response  amplitude  following  two 
opposite  trends  above  and  below  the  ambient  temperature:  decreasing  response  with  increasing 

temperature  between +20  °C  <  T  <+60°C,  and  increasing  response  with  increasing  temperature 


between -40  °C  <T  <  +20  °C .  The  reason  for  the  opposite  trends  resides  on  the  competing  roles 
of  the  piezoelectric  coefficients  and  the  dielectric  permittivity  terms,  as  noted  in  the  earlier  work 
on  isotropic  plates  by  Lanza  di  Scalea  and  Salamone  (2008). 


+a*-c 

♦so-c 


100  190  200  300 

flkMi] 


Figure  9.  MFC  pitch-catch  response  to  so  Lamb  mode  in  a  [0/±45/0]s  T700/5208  CFRP 
composite  laminate  under  changing  temperature  in  the  ranges  (a)-(b)  above  ambient  (+20°C  to 
+60°C),  and  (c)-(d)  below  ambient  (-40°C  to  +20°C).  Comparison  is  between  (a)-(c)  model  and 
(b)-(d)  experiment. 
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2.4.  Detection  Of  Bond  Condition  In  CFRP  Skin-To-Spar  Joint  Under  Changing 
Temperature 

This  section  discusses  the  detection  of  adhesive  bond  defects  in  a  specimen  simulating 
the  CFRP  skin-to-spar  joint  of  a  UAV.  The  test  joint  was  fabricated  by  bonding  a  [0/i  45/0]s 
T700/5208  CFRP  laminate  (dimensions  330  mm  x  330  mm  x  1.067  mm)  on  a  woven  T800/924 
CFRP  square  tube  (outer  dimensions  50.8  mm  x  50.8  mm,  wall  thickness  5.23  mm).  A  two-part 
Hysol  9394  epoxy  adhesive  was  used  as  the  bonding  agent.  Bonding  was  done  such  that  the  0- 
deg  fiber  direction  of  the  laminate  ran  along  the  lengthwise  direction  of  the  spar.  The  test  joint, 
showing  the  MFC  actuator/sensor  pairs,  is  shown  in  Figure  10. 


Figure  10.  CFRP  skin-to-spar  joint  specimen:  (a)  simulated  bond  defects  and  attached  MFC 
transducers;  (b)  in  environmental  chamber  for  temperature  tests. 

Defects  were  engineered  at  the  skin-spar  bondline  prior  to  joining.  The  majority  of  the 
bond  was  comprised  of  a  properly-mixed  epoxy,  representing  the  well-bonded  region.  A  poorly- 
cured  region  was  created  by  improperly  mixing  the  resin  and  hardener  resulting  in  an 
approximate  50%  degradation  in  longitudinal  and  shear  stiffness  compared  to  the  properly-cured 
case.  The  degradation  in  stiffness  was  verified  by  ultrasonic  pulse-echo  measurements  of 
longitudinal  and  shear  waves  conducted  on  epoxy  samples.  In  addition,  two  sizes  of  Teflon 
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release  film  inserts  (12.7  mm  x  12.7  mm  and  25.4  mm  x  25.4  mm)  were  added  representing 
disbonded  regions.  The  release  film  was  expected  to  severly  degrade  the  shear  stiffness  of  the 
bond  at  those  locations. 

It  was  previously  shown  (Lanza  di  Scalea  et  al,  2007)  that  the  strength  of  transmission  of 
the  So  mode  through  the  joint  is  sensitive  to  the  presence  of  these  bond  defects  at  room 
temperature.  That  study  also  showed  that  maximum  sensitivity  to  bond  state  was  achieved  at 
“mode  coupling”  conditions  (corresponding  to  a  large  transfer  of  ultrasonic  energy  through  the 
thickness  layers  of  the  joints).  Mode  coupling  occurred  for  So  at  frequencies  around  200  kHz  for 
the  subject  joint.  This  frequency  value  also  corresponded  to  an  efficient  wave  transduction 
through  wavelength  tuning  for  the  given  MFC  transducer  patches  on  the  [0/  — 45/0]s  CFRP  plate 
according  to  the  model  presented  earlier. 

The  present  work  adds  the  temperature  variation  to  the  study  by  Lanza  di  Scalea  et  al. 
(2007).  It  was  discussed  in  the  previous  sections  that  temperature  excursions  can  produce 
significant  changes  in  guided-wave  amplitude.  Hence  an  effective  on-board  monitoring  system 
must  be  able  to  discriminate  changes  due  to  bond  defects  from  those  due  to  temperature 
variations.  This  problem  was  here  addressed  as  statistical  pattern  recognition  by  a  multivariate 
Outlier  Analysis. 

2.4. 1 .  Experimental  Procedure 

A  total  of  six  MFC  pairs  were  installed  to  perform  pitch-catch  test  through  the  bonded 
region  of  the  test  joint  (Figure  10a).  Three  pairs  monitored  three  different  areas  of  the  properly- 
cured  bond  for  statistical  significance.  The  other  pairs  monitored  each  of  the  three  defected 
regions,  namely  the  poorly-cured  region  and  the  two  (large  and  small)  disbonded  regions.  The 
joint  was  tested  in  the  environmental  chamber  (Figure  10b)  between  -40  °C  and  ^60  ®C  in  10  ®C 
steps. 

Ultrasonic  signals  were  acquired  by  all  six  pairs  at  205  kHz  and  time  gated  to  isolate  the 
predominant  so  mode.  Eight  features  related  to  the  strength  of  transmission  of  the  wave  were 
computed  from  the  measurements.  These  features  were  the  Root-Mean-Square,  the  Variance,  the 
Kurtosis,  the  Peak-to-Peak  and  the  Peak  of  the  time-domain  signal;  the  Peak  and  the  Area  of  the 
signal’s  Fourier  transform  spectrum;  and  the  Area  of  the  signal’s  Hilbert  transform  spectrum. 


28 


These  features  were  assembled  into  an  8-dimensional  Feature  Vector  for  each  measurement 
point. 

A  multivariate  Outlier  Analysis  was  conducted  by  computing  the  known  Mahalanobis 
Squared  Distance  (MSD)  as  the  discordancy  metrics  (Worden  et  al.  2002): 


(20) 


where  x  is  the  current  Feature  Vector,  x  is  the  mean  vector  of  the  baseline,  S  is  the  covariance 
matrix  of  the  baseline,  and  T  represents  a  transposed  matrix.  If  the  MSD  value  exceeds  a 
predetermined  threshold,  x  is  an  outlier  i.e.  statistically  different  from  the  baseline  distribution. 
The  threshold  was  here  computed  by  a  Monte  Carlo  simulation  corresponding  to  a  99.9% 
confidence  interval,  i.e.  1  of  every  1000  observations  would  result  in  a  false-positive 
classification. 

Regarding  the  baseline  distribution,  two  different  cases  were  considered.  In  one  baseline, 
the  set  incorporated  the  measurements  taken  for  the  three  MFC  pairs  probing  the  properly-cured 
bond  over  the  entire  temperature  range  (-40  ®C  to  60  °C).  This  is  an  “ideal”  baseline,  because  it 
would  simply  require  storing  one  set  of  measurements  to  perform  statistical  defect  detection 
during  the  entire  range  of  aircraft  operations.  However,  this  is  also  an  “ambitious”  baseline 
because  it  limits  the  ability  of  the  pattern  recognition  to  properly  flag  the  defects.  Hence  a  second 
baseline  choice  was  also  considered,  whereas  the  full  temperature  range  was  partitioned  into  two 
ranges,  -40  °C  to  20  °C  (“below  ambient”)  and  20  °C  to  60  °C  (“above  ambient”).  In  this  second 
case  the  baseline  distribution  will  have  to  be  changed  depending  on  whether  the  aircraft  is  flying 
or  it  is  stored.  This  particular  choice  of  baseline  partition  was  consistent  with  the  two  different 
trends  observed  for  the  pitch-catch  response  in  the  composite  plate 

2.4.2.  Results 

Figure  1 1  shows  the  results  of  the  univariate  discordancy  test  based  on  only  one  signal 
feature.  Root  Mean  Square  is  shown  in  Figure  1  la  and  Peak-to-Peak  is  shown  in  Figure  1  lb.  The 
baseline  set  incorporates  measurements  taken  for  the  three  properly-cured  regions  over  the  entire 
temperature  range  (-40  °C  to  60  ®C).  The  99.9%  threshold  is  shown  as  a  horizontal  straight  line. 
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Each  point  corresponds  to  a  measurement  at  a  specific  temperature  and  joint  condition.  It  can  be 
seen  that  the  univariate  analysis  is  not  able  to  properly  detect  the  bond  defects  as  outliers,  with 
the  exception  of  a  few  measurements  of  the  poorly-cured  bond.  Hence  in  this  analysis 
temperature  effects  mask  the  effects  of  the  defects. 
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Figure  11.  Univariate  Outlier  Analysis:  results  of  discordancy  test  for  (a)  Root  Mean  Square 
signal  feature;  (b)  Peak-to-Peak  signal  feature 
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The  results  of  the  multivariate  analysis  which  uses  the  8-dimensional  Feature  Vector  are 
shown  in  Figure  12.  The  full  baseline  distribution  (-40  °C  to  60  °C)  is  again  used  in  this  figure.  It 
is  observed  that  the  discriminatory  power  of  the  multivariate  case  is  greatly  improved  relative  to 
the  univariate  case,  with  many  of  the  measurements  taken  at  the  defected  regions  being  properly 
classified  as  outliers  (MSD  above  threshold).  However,  some  data  points  of  the  damage  cases, 
particularly  the  “small”  disbond  case,  are  improperly  flagged  as  inliers  i.e.  they  represent  false 
negatives. 
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Figure  12.  Multivariate  Outlier  Analysis:  results  of  discordancy  test  with  baseline  distribution 
from  full  temperature  range  (-40  °C  to  60  ®C) 


A  further  improved  defect  detection  performance  is  obtained  when  the  baseline  is 
partitioned  between  the  “below  ambient”  temperatures  (Figure  13a)  and  “above  ambient” 
temperatures  (Figure  13b).  In  these  analyses  all  measurements  at  defected  regions  are  properly 
flagged  as  anomalies  (no  false  negatives).  Furthermore,  all  measurements  at  properly-cured 
regions  are  correctly  flagged  as  inliers  (no  false  positives). 
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Figure  13.  Multivariate  Outlier  Analysis:  results  of  discordancy  test  with  baseline  distribution 
partitioned  between  (a)  below  ambient  temperature  (-40  °C  to  20  °C),  and  (b)  above  ambient 
temperature  (20  °C  to  60  °C). 


These  results  suggest  that  it  is  possible  to  detect  the  bond  defects  considered  in  this  study 
within  the  temperature  excursion  of  normal  aircraft  operations  and  in  a  statistically  robust 
manner  once  (1)  guided-wave  measurements  are  taken  at  mode-frequency  points  that  are 
particularly  sensitive  to  the  bond  layer  (so  at  ~200  kHz  for  the  subject  joint),  (2)  several  features, 
rather  than  any  single  feature,  are  extracted  from  the  guided-wave  measurements  to  perform  a 
multivariate  Outlier  Analysis,  and  (3)  a  partition  of  the  baseline  distribution  is  performed  to 
improve  the  discriminatory  power  of  the  statistical  algorithm. 


2.5  Discussion  and  Conclusions 


The  first  part  of  the  section  discussed  a  model  to  predict  the  amplitude  spectrum  of  pitch- 
catch  guided-wave  signals  in  composite  plates  probed  by  MFC  piezocomposite  transducers.  The 
model  accounts  for  the  transducer  piezo-mechanical  properties,  the  transducer-plate  interaction 
through  shear-lag  effects,  and  the  substrate  plate  properties,  all  of  which  can  be  parametrized 
with  temperature.  Theoretical  and  experimental  results  obtained  on  a  [0/±45/0]s  CFRP  laminate 
between  -40  °C  and  60  °C  show  substantial  changes  in  the  predominant  so-mode  amplitude  (as 
high  as  40%  change  from  the  ambient-temperature  value).  Little  change,  instead,  is  seen  in  the 
wavelength  tuning  points  (maxima  and  zeros  in  the  response).  It  was  also  observed  that  the 
signal  amplitude  follows  two  opposite  trends  below  and  above  ambient  temperature. 
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The  second  part  of  the  section  has  applied  these  findings  to  the  detection  of  bond  defects 
in  a  CFRP  skin-to-spar  adhesive  joint.  The  defects  considered  were  a  region  of  poorly-cured 
adhesive  (50%  degradation  in  stiffnesses)  and  two  different  sizes  of  disbonds.  The  study  used  an 
So  ultrasonic  frequency  (~200  kHz)  which  corresponded  to  a  wavelength  tuning  peak  for  the 
MFC  transducers,  as  well  as  to  a  mode-coupling  condition  for  the  multilayered  joint  structure. 
The  strength  of  transmission  of  this  frequency  was  used  to  monitor  bond  state  in  the  (-40  °C  to 
60  °C)  temperature  range.  It  was  found  that  a  multivariate  Outlier  Analysis  using  8  representative 
signal  features  can  be  successful  in  detecting  the  bond  defects  within  the  temperature  excursions. 
Correct  defect  classification  was  obtained  when  the  baseline  distribution  for  the  analysis  was 
partitioned  in  two  regions,  below  ambient  temperature  (-40  °C  to  20  °C)  and  above  ambient 
temperature  (20  °C  to  60  °C).  In  an  on-board  system,  this  suggests  that  the  defects  examined  can 
be  successfully  detected  by  toggling  the  baseline  between  only  two  operational  conditions:  in¬ 
flight  and  storage. 

There  are  several  assumptions  in  the  present  work.  The  model  assumes  plane  waves,  thin 
transducers  and  neglects  damping  in  the  substrate  plate.  Also,  the  results  on  the  bond  defect 
detection  only  apply  to  the  specific  CFRP  skin-to-spar  joint  tested.  Different  choices  of  wave 
frequency.  Feature  Vector  dimension,  and  baseline  distribution  will  likely  be  required  for 
different  applications. 

3  GLOBAL-LOCAL  MODELING/SLMULATION  APPROACH  FOR  QUANTITATIVE 
DEFECT  ASSESSMENT  USING  GUIDED  WAVES 

3.1.  Theory 

The  geometry  of  the  problem  at  hand  is  depicted  in  Figure  14  where  the  material  inside  the 
bounded  surface  has  properties  different  from  the  outside  due  to  the  presence  of  local  irregularities 
as  a  defect.  In  order  to  best  illustrate  the  hybrid  method,  assume  that  a  plane  time-harmonic  wave  is 
incident  on  the  FEM  boundary  in  Figure  14.  The  region  inside  the  boundary  is  meshed  and  analyzed 
using  standard  finite  element  analysis;  the  region  outside  the  boundary  is  represented  by  the  linear 
combination  of  a  set  of  global  wave  functions,  here  calculated  from  the  SAFE  method.  The 
continuity  of  displacements  and  tractions  on  the  mesh  boundary  must  be  satisfied  through  proper 
choice  of  the  amplitudes  of  the  wave  functions. 
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Figure  14.  Geometry  of  a  general  wave  scattering  problem  by  Global- Local  method. 

3.1.1.  Semi-Analytical  Finite  Element  Method 

The  SAFE  method  is  employed  to  calculate  the  guided  wave  modeshapes  O  and  eigenvalues 
^  for  any  infinite ly-uniform  cross-section.  These  modeshapes  and  eigenvalues  are  used  to  calculate 
the  forces  and  displacements  on  the  boundary  separating  the  Global  from  the  Local  region.  In  the 
SAFE  method,  at  each  frequency  co  a  discrete  number  of  guided  modes  is  obtained.  For  the  given 
frequency,  each  mode  is  characterized  by  a  wavenumber,  and  by  a  displacement  distribution  over 
the  cross-section.  For  a  Cartesian  reference  system,  the  waveguide  cross-section  is  set  in  the  y-z' 
plane  while  the  x-axis  is  parallel  to  the  waveguide  length  (see  Figure  15). 

i  Mono-dimensional 
Element 


Figure  15.  Plate  waveguide  cross-section  in  Cartesian  coordinates  for  SAPE. 
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Subdividing  the  cross-section  via  finite  elements,  the  approximated  displacement  at  a  point,  m,  is 
given  as; 


(21) 


where  N(y,z)  is  the  matrix  of  the  shape  functions,  U"  is  the  nodal  displacement  vector  for  the  e'^ 
element,  t  is  the  time  variable  and  i  the  imaginary  unit.  It  can  be  noted  that  the  displacement  is 
described  by  the  product  of  an  approximated  finite  element  field  over  the  waveguide  cross-section 

with  exact  time  harmonic  fimctions,  in  the  propagation  direction,  ‘x’.  The  compatibility 

and  constitutive  equations  can  be  written  in  synthetic  matrix  forms  as: 

£  =  Du,  a  =  C'e  (22) 

where  f  and  dare  the  strain  and  stress  vectors  respectively,  D  is  the  compatibility  operator,  and  C” 
is  the  stiffness  tensor.  More  details  on  the  compatibility  operator  can  be  found  in  Hayashi  et  al.  [8] 
and  Bartoli  et  al.  [9].  The  principle  of  virtual  work  with  the  compatibility  and  constitutive  laws. 
Equations  (21)  and  (22),  leads  to  the  following  energy  balance  equation 

J  du'^tdT  =  J  5u'^[pu)dV  +  J  5[uD'^  C' DudV  (23) 

where  F  is  the  waveguide  cross-section  area,  V  is  the  waveguide  volume,  t  is  the  external  traction 

vector  and  the  overdot  implies  time  derivative.  The  finite  element  procedure  reduces  Equation  (23) 
to  the  set  of  algebraic  equations; 


(24) 


where  the  subscript  2M  indicates  the  dimension  of  the  problem  with  M,  the  number  of  total  degrees 
of  freedom  of  the  cross-sectional  mesh.  Details  on  the  complex  matrices  [A]  and  [B]  can  be  found 
in  Bartoli  et  al.  (2006).  By  setting  p=Q  in  Equation  (24),  the  associated  eigenvalue  problem  can  be 
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solved  as  ^((o).  For  each  frequency  co,  2M  complex  eigenvalues  ^  and  2M  complex  eigenvectors 
<I>m  are  obtained.  The  solution  is  symmetric,  i.e.  for  each  pair  representing  a  forward 

guided  mode,  a  pair  exists  representing  the  corresponding  backforward  mode  (Hayashi  et  al.  [8]). 
The  first  M  components  of  Om  describe  the  cross-sectional  mode  shapes  of  the  w-th  mode.  Once  4m 
is  known,  the  dispersion  curves  can  be  easily  computed.  The  phase  velocity  can  be  evaluated  by  the 
expression  Cph=(olkreai,  where  kreai  is  the  real  part  of  the  wavenumber.  The  imaginary  part  of  the 
wavenumber  is  the  attenuation,  att=  k^ag,  in  Nepers  per  meter. 


3.1.2  Global-Local  Method 

The  global  functions,  represented  by  On,  are  the  displacements  that  satisfy  the  governing 
equations  in  the  absence  of  the  discontinuity.  The  displacements  at  the  mesh  boundary  are  written 
as  a  sum  of  the  contribution  due  to  the  incident  wave  and  that  due  to  the  scattered  waves.  This 
displacement  is  then  further  used  to  satisfy  the  continuity  conditions  at  the  boundary.  The 
displacements  at  the  left  boundary  are  written  as: 


u=ii  +u  .  =  A.  O*  e'^  ^ +  'y  A  .O  e 

I  incident  scattered  in  in  j  j 


y=i 


I  in  in 


(25a) 

(25b) 


where  the  time-harmonic  dependence  is  implicit. 

Rattanawangcharoen  et  al.  [4]  give  a  detailed  description  of  matrices  [G"]  and  [D'].  In 
Equation  (25),  represents  the  modeshape  of  the  incident  waves  travelling  to  the  right,  O" 
represents  the  modeshapes  of  the  reflected  waves  travelling  to  the  left;  and  represent 

the  amplitude,  modeshape  and  wavenumber  of  the  incident  wave  respectively;  AJ  represents  the 
unknown  amplitude  of  the  reflected  waves;  superscript  represents  a  wave  travelling  in  the  right 

and  a  superscript  represents  a  wave  travelling  in  the  left  direction;  represents  the  absolute 
value  of  the  distance  of  the  left  boundary  from  the  center  of  the  meshed  region;  N  is  the  total 
number  of  modes  considered.  In  the  same  way,  the  displacements  at  the  right  boundary  are  written 
as: 
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[G‘][0-] 


«,  =  «  w  ,  +  "  „  ,  =  A  '  +  y  A^(i>^e 

I  incident  scattered  in  in  j  j 

7=1 


(26a) 

(26b) 


In  Equation  (26),  represents  the  unknown  amplitude  of  the  transmitted  waves  and  represents 

the  absolute  value  of  the  distance  of  the  right  boundary  from  the  center  of  the  meshed  region. 

The  forces  at  the  boundaries  can  be  calculated  by  computing  the  consistent  nodal  force 
vectors  due  to  displacements  (5b)  and  (6b).  The  nodal  forces  at  the  left  boundary  are: 

(27) 


Fj  is  the  consistent  nodal  force  vector  calculated  from  the  displacement  vector  d>j.  The  nodal  forces 
at  the  right  boundary  are: 


fr=  *[F'][D-] 


(28) 


3.1.3  Finite  Element  Model  for  Local  Region 

The  conventional  discretization  process  in  the  finite  element  method  yields: 

(29a) 
(29b) 

(29c) 


Here,  [/f;],  and  [W/]  are  global  stiffness  and  mass  matrices  respectively;  co  is  the  circular 
frequency;  {q;}  and  {q^}  are  the  nodal  displacement  vectors  corresponding  to  interior  and 
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boundary  nodes  respectively;  is  the  interaction  force  vector  at  the  boundary  nodes;  and  q” 
denotes  complex  conjugate  of  q . 

Forces  on  the  boundary  nodes  due  to  the  entire  FE  region,  when  displacement  is  imposed 
only  on  the  boundaries,  can  be  calculated  by  performing  static  condensation  on  [5]  such  that; 


(30) 


3.1.4  Global  Solution 

The  nodal  forces  on  the  boundary  due  to  displacements  (5b)  and  (5b),  as  contributed  from 
the  FE  portion,  are: 


hEM 

houndarv 


(31) 


For  traction  continuity,  from  Equations  (25),  (26),  and  (21): 


i 


=  h 


FEM 

boundary 


(32) 


Rearranging  the  terms,  the  following  final  equation  is  obtained: 

([s»][c]-[f]){o}={/::}-[^;,]{9'} 


' - 

m  tn  1 

A  0U^  '"  ’ 

m  tn 

(33a) 

(33b) 


(33c) 
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(33d) 


{D}% 

=(  {^-r  {^-i 

[^]  = 

[F-]  0 

0  [r] 

(33e) 


0  [O*] 


Equation  (33a)  is  solved  using  least  squares  to  calculate  the  coefficients  Aj  ,  and  A* .  Reflection 
coefficients  R.j  and  transmission  coefficients  (due  to  the  /'*  incident  mode)  can  be  found  by: 


T.= 

ij 


r,  =  1  + 

IJ 


(35) 


3.1.5  Energy  Conservation 

Energies  are  carried  only  by  the  propagating  modes,  which  include  incident,  reflected  and 
transmitted  modes.  The  time  averaged  values  of  energy  flux  associated  with  the  reflected  and 
transmitted  modes  through  the  cross-section  is  given  by: 


r.=03 
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where 
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(37) 


h 

a  =  \  F^q  zdz 

j  J  j 
0 


is  the  participation  factor  for  the  j*  propagating  mode. 


3.1.6.  Forced  Solution  from  SAFE  and  Global-Local 

SAFE  can  be  used  to  compute  the  displacement  time  history  at  any  point  due  to  an  arbitrary 
loading  over  a  cross-section.  The  displacement  field  is  written  as  a  linear  combination  of  the  wave- 
modes  present  in  the  structure  (Hayashi  et  al.  [8]).  Mathematically, 

(38) 

m=l 

where,  0^“'"  is  the  upper  half  of  the  right  propagating  eigenvector,  is  the  corresponding 
eigenvalue,  x  is  the  coordinate  location  of  the  point  in  consideration  and  is  the  coordinate 
location  of  the  applied  force;  represents  the  individual  contribution  of  the  mode. 

The  displacement  vector,  U,  is  calculated  from  Equation  (38)  for  discrete  frequencies  and  an 
inverse  Fourier  Transform  is  used  to  invert  the  frequency  domain  to  the  time  domain.  For 
calculating  time  histories  from  the  Global-Local  method,  the  time  domain  forcing  function  is  first 
converted  to  its  corresponding  frequency  domain  representation.  For  each  frequency  component. 
Equation  (38)  is  used  to  express  the  displacement  field  due  to  the  force  as  a  linear  combination  of 
the  displacement  eigenvectors,  .  For  each  such  eigenmode  contribution,  Global-Local  is  then 

used  to  calculate  the  scattered  field.  This  scattered  field  contributes  to  the  total  displacement  at  the 
point  of  interest.  The  final  time  history  is  found  by  inverting  the  sum  of  all  such  contributions  for  all 
the  different  excitation  frequencies  under  consideration. 
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3.2.  Results 


3.2. 1  ■  Reflection  for  Notched  Aluminum  Plate 

The  first  case  examined  is  a  semi-infinite  isotropic  plate  with  a  notch  defect.  The  plate 
modeled  was  aluminum,  1 .58  mm  (1/16  in)  in  thickness.  The  first  notch  considered  was  0.16  mm  in 
width  and  0.32  mm  in  depth  (Figure  16a).  A  finite  region  surrounding  the  notch  (local  domain)  and 
extending  into  the  thickness  of  the  plate,  was  discretized  by  rectangular,  four-node  finite  elements. 
A  total  of  ten  elements  were  used  in  the  thickness  of  the  plate,  and  the  mesh  was  kept  uniform 
throughout  the  local  domain.  SAFE  was  used  for  the  global  domain.  The  zero-order  symmetric 
mode,  Sq  ,  was  considered  as  the  incident  mode. 
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Figure  16.  Results  of  So  wave  scattering  from  defects  in  a  plate  from  Global-Local  models  in 
reflection,  (a)  Energy  spectra  of  reflected  So  and  Ao  (mode  conversion)  from  sharp  notch;  (b) 
Energy  spectra  of  reflected  So  and  Ao  (mode  conversion)  from  shallow  notch. 


The  results  are  shown  in  the  plot  of  Figure  16a  in  terms  of  reflected  energy  relative  to  the 
incident,  unit  Sq  energy.  As  expected,  the  notch  causes  Sq-A^  mode  conversion  since  it  breaks  the 
symmetry  of  the  waveguide.  As  a  check  for  the  accuracy  of  the  analysis,  energy  conservation  is 
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preserved  as  the  reflected  Sq  and  energies  sum  to  the  incident  Sq  energy.  Mode  conversion  is 

found  more  pronounced  at  particular  frequency  values,  namely  500  kHz,  730  kHz,  880  kHz  and  950 
kHz. 

A  more  shallow  discontinuity  was  examined  next.  This  was  a  notch  extending  1.6  mm  in 
width  and  0.16  mm  in  depth,  which  more  closely  represented  corrosion  damage  (Figure  16b).  It  can 
be  seen  that  the  resulting  reflection  spectra  are  different  for  the  two  defect  cases.  For  example,  the 
strong  Sg-A^  reflection  at  730  KHz  in  the  sharp  defect  is  now  completely  absent  in  the  shallow 

defect.  Maximum  Sg-A^  reflection  for  the  shallow  defect  is,  instead,  found  at  550  KHz  where  zero 

reflection  is  seen  in  the  sharp  defect.  This  suggests  the  potential  for  quantitative,  rather  than 
qualitative  defect  detection. 


3.2.2.  Reflection  and  Transmission  for  Notched  Aluminum  Plate 

In  this  case  the  global-local  method  was  used  to  predict  the  scattered  field  in  both  reflection 
and  transmission.  The  plate  modeled  was  aluminum,  1.6  mm  in  thickness  with  a  half  thickness 
notch  in  the  middle.  Transmitted  and  reflected  Sg  and  Ag  modes  were  computed  for  an  incident  Sg 
mode.  Energy  conservation  was  again  checked  to  test  the  validity  of  the  method  (Figure  17). 


SAFE  Global  2-D  K  local _ SAFE  Global 

So  Incident* -  7  mm  — h  /2*0^mmhofe  so  transmitted 

- - - 1  ^  I - n - 


Soreflectedi 
Ao  reflected 


16  mm 


1.S  mm 

t 


Ao  transmitted 


Frequency  kHz 


Figure  17.  Results  of  So  wave  scattering  in  a  plate  from  Global-Local  model  in  reflection  and 
transmission. 
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3.2.3.  Reflection  and  Transmission  for  Layered  Composite 

This  example  examined  reflection  and  transmission  in  the  more  complex  case  of  a  layered 
composite  plate.  The  plate  modeled  was  1mm  thick  and  consisted  of  8  layers  of  unidirectional 
carbon-epoxy  laminae  in  a  stacking  sequence  of  [0/+45/-45/0]s.  Material  damping  was  ignored  in 
the  present  study.  A  defect  of  length  3mm,  consisting  of  the  removal  of  certain  layers,  was 
introduced  at  various  locations  across  the  cross-section  of  the  laminate  (Figure  18).  The  material 
density  of  each  lamina  was  p=1530  kg/m^  and  the  elastic  properties  in  the  principal  direction  of 
material  symmetry  are  given  in  Table  4,  where  1  is  the  fiber  direction,  2  is  the  direction 
perpendicular  to  the  fibers  in  the  laminate  plane,  and  3  is  the  through-thickness  direction.  For  waves 
propagating  along  a  direction  x  oriented  at  any  angle  0  with  respect  to  the  fiber  direction,  1,  the 
SAFE  model  simply  requires  the  rotation  of  the  stiflfiiess  matrix  of  each  lamina  through 


(39) 


where  C  is  the  stiffness  matrix  in  the  lamina’s  principal  direction,  and  R/  and  R2  are  the  following 
rotation  matrices: 
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with  m  =  cos(0)  and  n  =  sin  (0) .  The  governing  eigenvalue  problem  of  SATE  in  Equation  (24)  is 
then  solved  by  using  the  rotated  stiffness  matrix  in  the  constitutive  relations. 


c„ 

135 

C,2 

5.7 

C,3 

5.7 

C22 

14.2 

C23 

8.51 

C33 

14.2 

C44 

2.87 

C55 

4.55 

C66 

4.55 

Table  4:  Elastic  properties  (GPa)  for  the  CFRP  lamina. 

Figure  18  shows  the  scattered  fields  for  an  Sq  incidence  over  the  defected  layered 

composites.  Three  different  defects  are  considered.  Scattered  fields  from  the  removal  of  the  top  2 
layers  (Figure  18a)  are  compared  to  the  removal  of  the  middle  two  layers  (Figure  18b),  and  the 
removal  of  the  top  four  layers  (Figure  18c). 

It  can  be  seen  that  the  reflection  strength  increases  as  the  defect  is  moved  from  the  top  of  the 
plate  (Figure  18a)  to  the  center  of  the  plate  (Figure  18b).  This  can  be  explained  by  the  energy 
distribution  in  the  fundamental  symmetric  mode  over  the  cross-section  of  the  plate.  Figure  19  plots 
the  Poynting  vector  for  the  mode  at  the  frequency  of  ~400kHz  across  the  cross-section  of  the 
composite  plate,  where  the  Poynting  vector  is: 


P  =  - 


\_ 

2 


Re  {uu  ^ 


(41) 


where  a  is  the  classical  stress  tensor,  and  u*  is  the  complex  conjugate  of  the  particle  velocity 
vector. 
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Figure  18.  Results  of  So  wave  scattering  from  defects  in  a  [0/-^45/-45/0]s  CFRP  laminate  in 
reflection  and  transmission,  (a)  removal  of  top  two  layers;  (b)  removal  of  central  two  layers;  (c) 
removal  of  top  four  layers. 


^^3  Magnitude  of  Poynting  Vector  (a.u.) 


Figure  19.  Poynting  vector  for  So  at  400  KHz  in  the  [0/+45/“45/0]s  CFRP  laminate. 
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Figure  19  shows  that  the  energy  carried  by  the  Sq  mode  in  the  direction  of  propagation  is 

concentrated  more  around  the  center  of  the  cross-section  of  the  laminate  than  at  its  surface. 
Therefore,  it  is  expected  that  this  mode  would  be  more  sensitive  to  a  central  defect  compared  to  a 
defect  at  the  top  or  the  bottom  of  the  laminate.  Figure  18a  and  18b  accordingly  show  an  increased 
reflection  of  the  mode  when  the  middle  2  layers  are  removed  as  compared  to  the  case  when  the 

top  2  layers  are  removed.  In  addition,  mode  conversion  is  completely  absent  in  Figure  1 8b  given  the 
symmetry  of  the  defect  with  respect  to  the  cross-section;  alternatively,  considerable  mode 
conversion  of  Sg  into  is  seen  in  Figure  1 8c  because  the  defect  breaks  the  symmetry  of  the  cross- 

section.  As  expected,  Global-Local  predicts  a  significant  increase  in  scattering  as  the  defect  size 
increases  from  2  layer  thicknesses  to  4  layer  thicknesses  (Figure  18c). 

Ao  incidence  was  studied  next  for  the  same  composite  plate.  Figure  20  shows  the  scattered 
fields  for  the  same.  As  in  the  previous  problem,  three  different  defects  are  considered  (Figure  1 8  a- 
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Figure  20.  Results  of  Ao  wave  scattering  from  defects  in  a  [0/+45/-45/0]s  CFRP  laminate  in 


reflection  and  transmission,  (a)  removal  of  top  two  layers;  (b)  removal  of  central  two  layers;  (c) 


removal  of  top  four  layers. 
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Figure  21  shows  the  magnitude  of  the  Poynting  vector  for  the  mode  at  540  kHz.  Similar 

to  the  case  of  Sg,  a  very  large  part  of  the  energy  is  concentrated  in  the  central  portion  of  the 

laminate’s  cross-section.  This  explains  the  strong  scattering  from  the  central  defect  (Figure  20b) 
compared  to  the  top  defect  of  the  same  size  (Figure  20a). 


Figure  21.  Poynting  vector  for  Ao  at  540  KHz  in  the  [0/-*-45/-45/0]s  CFRP  laminate. 

It  is  likely  that  the  finite  element  discretization  was  not  adequate  for  solving  the  Global- 
Local  problem  at  higher  frequencies.  This  might  explain  why  energy  conservation  for  both  Sg  and 

y4g  incidences  (Figs.  18  and  10)  is  not  completely  satisfied  at  these  frequencies.  Since  /Ig  in  general 
has  a  shorter  wavelength  compared  to  Sg ,  it  is  more  affected  by  the  finite  element  discretization. 
Moreover,  the  laminate  being  anisotropic,  there  is  a  coupling  of  the  shear  horizontal  mode  with  both 
the  Sg  and  Ag  mode.  This  allows  energy  from  the  latter  two  modes  to  escape  into  the  SN  modes. 
Despite  these  limitations,  the  underlying  physical  behavior  of  the  system  is  captured  well. 

3.2.4.  Experimental  Verification 

The  Global-Local  method  was  also  used  to  predict  the  time-domain  waveforms  in  the  CFRP 
[0/+45/-45/0]s  composite  plate  subjected  to  broadband,  pulsed  laser  excitation  (ND:YAG,  10  nsec. 
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pulse  duration).  The  response  was  experimentally  measured  for  the  case  of  a  plate  1.0  mm  in 
thickness  with  a  sensor-excitation  distance  of  20  cm.  The  defect  modeled  was  a  delamination  type 
defect  extending  in  the  thickness  direction,  in  the  middle  of  the  8-layered  composite  laminate.  Figs. 
22(a)  and  22(b)  show  the  comparison  of  the  theoretical  predictions  with  experimental  results.  The 
signals  are  normalized  to  their  respective  maximum  values.  It  can  be  noticed  that  both  the 
experimental  and  the  numerical  results  match  well  in  the  time  of  arrival  and  phase  characteristics  of 
the  signals.  Moreover,  a  decrease  in  the  damaged  signal  compared  to  the  undamaged  case,  as 
predicted  by  the  numerical  model,  is  confirmed  in  the  experimental  results.  Another  interesting 
result  is  that  both  predictions  and  experiment  show  a  velocity  change  between  the  undamaged  and 
damaged  cases.  This  is  likely  due  to  mode  conversions. 

3.3.  Discussion  And  Conclusions 

This  section  demonstrates  the  use  of  SAFE  in  a  Global-Local  framework  that  allows 
extracting  information  on  guided  waves  scattered  by  defects  in  isotfopic  and  layered  anisotropic 
plates  for  quantitative  structural  health  monitoring.  The  advantage  of  SAFE  is  its  ability  to  handle 
the  complex  geometry  of  layered  composites  in  a  computationally  efficient  manner.  The  method 
was  used  in  cases  of  notch-like  defects  in  aluminum  plates  and  delamination-like  defects  in  layered 
composites.  Both  reflected  and  transmitted  spectra  of  propagating  modes  were  predicted.  It  should 
be  noted  that  two  different  discretization  are  done  in  the  solution  of  the  problem.  These  affect  the 
accuracy  of  the  solution  in  different  limiting  regimes.  For  example,  the  thickness  profile  of  the  out 
of  plane  displacement  component  of  the  *5*0  mode  in  an  Aluminum  plate  varies  from  being 
approximately  linear  at  an  /c^  value  of  0.5  to  a  half  wavelength  sinusoid  at  f  d=  2.0  (Rose  1999). 
This  indicates  that  the  thickness  discretization  of  the  Global  part  in  Global-Local  should  be 
increased  to  adequately  model  the  increasing  nonlinear  thickness  variation  of  the  modeshape  as  the 
frequency  is  increased.  Higher  order  modes  exhibit  a  more  nonlinear  modeshape  compared  to  lower 
order  modes.  ,  for  instance,  has  more  wavelengths  in  its  displacement  thickness  profile  f  d=  5 

than  has  at  the  same  /  d  value.  Therefore,  in  general,  a  higher  order  mode  needs  a  finer  thickness 

discretization  at  the  same  frequency  than  a  lower  order  mode  in  order  to  better  approximate  the  true 
modeshape. 


48 


1*5 


Of 

-a 

E 

< 

T5 


1 

0*5 

0 

-0.5 

-1 


-Undainaged 
^  Damaged 


OOOi+00 

(a) 


5.00E-05 


l.DOE^^ 
Time  (s) 


l.SOf-Ofl 


2.00E-04 


(b)  T1"«(S) 


Figure  22.  Time  histories  (normalized  to  max.)  of  response  of  damaged  and  undamaged  [0/+45/- 
45/0]s  CFRP  laminate  subjected  to  broadband  laser  excitation,  (a)  Global-Local  predictions;  (b) 
Experimental  measurements. 


The  local  discretization  for  the  Finite  Element  part  of  Global  local  should  also  vary  with 
frequency  and  the  modes  under  consideration.  At  any  given  frequency,  modes  with  shorter 
wavelengths  require  a  finer  FE  discretization  for  good  approximation.  Modes  with  low  phase 
velocities  in  the  dispersion  curve,  therefore,  require  more  finite  elements  per  unit  propagation 
distance  for  an  accurate  representation  in  the  summation  solution.  For  isotropic  structures  like  an 
aluminum  plate,  since  all  the  modes  asymptotically  approach  a  constant  phase  velocity  (  Sq  and 
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converge  to  the  Rayleigh  surface  wave  velocity;  higher  order  modes  converge  to  shear  wave 
velocity),  the  FE  discretization  should  be  made  finer  as  higher  frequencies  are  considered. 

It  was  also  shown  that  the  wave  time  histories  can  be  extracted  from  these  models  for  an  arbitrary 
excitation.  The  results  shown  are  limited  to  a  2-D  study.  The  extension  of  the  work  to  3-D  cases  and 
to  adhesive  bond  assessment  should  be  the  subject  of  future  work. 

4.  A  TWO-STEP  DATA-DRIVEN  BOND  DEFECT  DETECTION  AND  CLASSIFICATION 
APPROACH 

The  research  completed  in  the  previous  sections  employed  fundamental  physics-based 
modeling  to  drive  algorithm  development.  One  of  the  limitations  of  the  subsequent  approach  is 
the  reliance  on  fairly  simple  structural  waveguide  geometries  except  through  the  cross-section. 
In  order  to  fiilly  augment  this  physics-based  approach,  this  section  employed  a  class  of  statistical 
pattern  recognition  algorithms  originally  developed  in  the  vibration  (lower-frequency  acoustic) 
domain  that  are  based  on  state  space  analyses  using  auto-  or  cross-prediction  methods.  These 
methods  used  chaotic  excitations  and  attractor-based  prediction  error  algorithms  and  have 
demonstrated  the  capacity  to  detect  damage  in  various  test  bed  structures  with  enhanced 
sensitivity  over  traditional  vibration-domain  analyses  (Nichols  et  al.  2003;  Todd,  Erickson  et  al. 
2004).  Unfortunately,  these  chaotic  excitations  are  also  low  frequency  in  nature  and  are  therefore 
unable  to  localize  damage  within  a  structure.  The  work  in  this  section  uses  the  benefits  of 
combining  ultrasonic  guided  waves  (small  length  and  time  scales)  with  chaotic  excitations 
(which  enable  and  enhance  pattern  recognition  techniques)  for  damage  detection  of  bolted  and 
adhesively-bonded  connections  via  bulk  insonification.  This  synthesis  of  techniques  is 
accomplished  by  shifting  the  energy  of  a  low-frequency  chaotic  process,  such  as  the  common 
Lorenz  signal,  into  the  ultrasonic  frequency  range  (>20  kHz)  and  launching  it  into  the  structure 
as  a  guided  wave.  The  chaotic  ultrasonic  signals  used  here  are  fundamentally  created  via  amplitude 
modulation,  i.e.,  by  multiplying  a  single  ultrasonic  frequency  tone  by  an  amplitude  envelope  that  is 
created  by  a  chaotic  process.  The  waveform  appears  as  a  narrowband,  chaotically-modulated  signal 
centered  at  the  same  central  frequency  as  the  original  ultrasonic  tone.  A  chaotic  waveform  is  able  to 
enhance  prediction  error  based  damage  detection  features  because  the  signal  is  deterministic,  and 
because  of  the  phase  space  diversity  associated  with  the  chaos.  The  deterministic  aspect  of  chaotic 
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inputs  was  used  in  the  context  of  SHM  by  Todd  et  al.  when  they  showed  that  an  attractor-based 
detection  feature  had  superior  performance  to  modal-based  features  at  detecting  loss  of  stiffness  in  a 
computational  spring-mass  system  when  excited  by  a  chaotic  waveform  (Todd  et  al.  2001).  A 
number  of  studies  involving  chaotic  excitation  were  then  completed  in  the  following  years  (Trickey 
et  al.  2002;  Nichols  et  al.  2003;  Nichols,  Todd,  and  Wait  2003;  Moniz  et  al.  2004;  Nichols  et  al. 
2004). 

4.1  Chaotic  Ultrasonic  Signal  Creation 

A  chaotic  ultrasonic  waveform  yn  is  created  by  simple  amplitude  modulation  of  an  ultrasonic 
carrier  as  detailed  in  the  following  Equation  (41): 

y„  =  sm{2KnfJ  +  (41) 

where  yn  is  the  modulated  waveform,  fs  is  the  sampling  frequency  (the  hardware  used  in  this 
dissertation  samples  maximally  at  4  MHz  for  a  high-fidelity  waveform),  fc  is  the  desired  carrier 
frequency  (typically,  in  the  hundreds  of  kHz  in  the  applications  in  this  dissertation),  d  is  the 
modulation  depth,  and  Xn  is  a  chaotic  waveform  produced  by  a  separate  chaotic  process.  Thus,  the 
waveform  Xn  appears  as  a  narrowband,  chaotically-modulated  signal  centered  at  the  desired  central 
carrier  frequency,  or  a  chaotic  ultrasonic  wave  (CU W).  In  this  work,  the  chaotic  process  was 
obtained  as  the  JC/  output  of  the  well-known  Lorenz  system 

•^1=  10(^2  -^l) 

X2=  +  28x,  -  ^2 )  (42) 

^3  =(^1^2  /3) 


There  is  nothing  unique  about  the  Lorenz  system  for  generating  chaotic  output;  any  system  capable 
of  producing  a  chaotic  output  is  suitable,  but  the  Lorenz  system  shown  in  Eq.  (42)  has  a  robust 
parameter  region  for  producing  chaotic  output  and  was  selected  for  this  study.  The  Lorenz  system  is 
integrated  using  a  time-step  R/fs,  where  /?  is  a  frequency  ratio  that  can  be  modified  to  change  the 
fundamental  time  scales  of  the  chaotic  signal.  This  dissertation  uses  several  values  of  R  that  affect 
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the  frequency  regime  in  which  the  power  of  the  chaotic  signal  lies.  This  chaotic  signal  is  normalized 
through  division  by  the  maximum  of  the  absolute  value  of  the  signal  so  that  the  values  range  from  - 
1  to  1.  The  modulation  depth  d  controls  signal  bandwidth,  and  if  it  is  greater  than  one,  the  resulting 
signal  will  be  over-modulated  and  will  result  in  a  phase  inversion  at  the  points  where  |c/*x„|>l. 
These  phase  inversions  would  be  detrimental  to  any  prediction  algorithm,  and  d  is  therefore 
restricted  to  the  range  0<r/<l.  Figure  23  shows  the  effect  of  changing  the  frequency  ratio  R  and  the 
modulation  depth  d  on  the  time  series  and  power  spectral  density  of  the  interrogating  wave. 
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Figure  23.  Modulated  chaotic  signal  power  spectral  density  and  time  series  using 
R  =  0.05  and  d  =  0.4  (upper  left),  R  =  0.05  and  c/  =  1 .0  (mid  left)  and  R  =  0.30  and  d  = 
1 .0  (upper  right). 


4.1.1  Dimensionality  of  Chaotic  Ultrasonic  Signals 

Recently,  methods  based  on  state-space  analysis  (of  reconstructed  dynamic  responses)  have 
been  used  as  a  means  of  change  detection  in  SHM  applications.  In  short,  geometric  objects  that  are 
akin  to  the  true  state-space  attractors  of  the  chaotic  waveform  are  reconstructed  by  embedding  time- 
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delayed  copies  of  the  original  time  series  into  higher  dimensions  to  create  a  pseudo-attractor  in 
state-space  from  which  various  correlation  or  prediction  metrics  may  be  computed  to  facilitate  the 
feature  extraction  part  of  the  pattern  recognition  paradigm.  Features  that  have  been  developed  for 
state- space  analysis  in  general  include  correlation  dimension  (Grassberger  and  Procaccia  1983; 
Logan  and  Mathew  1996),  dynamic  continuity  (Pecora  et  al.  1997),  nonlinear  prediction  error 
(Schreiber  1997),  generalized  interdependence  (Amhold  et  al.  1999),  time  scale  separation 
(Cusumano  and  Chatteijee  2000),  Lyapunov  exponents  (Trendafilova  and  Van  Brussel  2001),  local 
attractor  variance  ratio  (Todd  et  al.  2001),  parameter  drift  (Chateijee  et  al.  2002),  attractor 
dimension  (Nichols  et  al.  2003)  and  phase  space  warping  (Chelidze  and  Cusumano  2004).  The  idea 
is  that  these  features  may  be  compared  between  baseline  and  test  cases  to  determine  how  the  system 
response  has  evolved  due  to  damage. 

A  critical  issue  when  employing  such  state-space  embedology  methods  is  the  embedding 
dimension  choice  that  will  be  used  in  the  pseudo-attractor  reconstruction  (Sauer  et  al.  1991).  Every 
orbit  in  a  multivariate  space  has  a  characteristic  fractional  dimension  that  describes  the  way  in 
which  the  number  of  points  within  a  sphere  of  radius  r  scales  as  the  radius  shrinks  to  zero.  There  are 
various  methods  for  determining  the  fractional  dimension  of  a  signal  (Abarbanel  1996).  The 
reconstruction  theorem  states  that  at  a  minimum  the  embedding  dimension  must  be  at  least  twice  the 
fractional  dimension  of  a  signal  to  ensure  that  the  signal  is  completely  unfolded  in  state-space  and 
will  behave  dynamically  similar  to  the  true  system.  However,  this  is  an  upper  bound  to  the  choice  of 
embedding  dimension,  and  it  is  usually  possible  to  embed  in  a  smaller  dimension  than  this  sufficient 
criterion  dictates. 

In  SHM  applications,  it  is  desirable  to  keep  the  embedding  dimension  as  low  as  possible  to 
maintain  damage  detection  resolution  (statistically,  many  observations  of  low-dimensional  features 
are  superior  to  fewer  observations  of  high-dimension  features  for  classification  problems)  and  to 
reduce  computational  demands,  where  appropriate.  A  low  embedding  dimension  is  especially 
necessary  if  these  attractor-based  methods  are  to  be  employed  in  real-time  monitoring  systems, 
given  the  processing  efforts  required  to  execute  them,  e.g.  nearest  neighbors  searching.  In  active 
sensing,  where  the  interrogating  ultrasonic  waveform  is  prescribed  with  a  certain  dimension 
associated  with  the  chaotic  process  used  to  create  it,  the  signal  is  processed  and  received  through 
various  data  acquisition  and  filtering  components  that  themselves  can  be  considered  dynamic 
systems.  It  is  known  that  these  actions  can  affect  the  dimension  of  a  process  (Pecora  and  Carroll 
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1996).  Therefore,  a  dimensionality  study  is  undertaken  to  examine  the  possible  effects  that  a 
standard  DAQ  process  has  on  signal  dimension  through  the  launching  process.  The  expected 
filtering  due  to  waveform  interaction  with  the  structure  itself  may  very  well  cause  dimension 
increase  as  well,  but  undesired  additional  dimensional  increase  due  to  the  DAQ,  conditioning,  or 
signal  generation  steps  must  be  minimized  for  optimal  efficacy  of  the  overall  approach. 

A  standard  data  acquisition  (DAQ)  process  can  be  represented  using  the  block  diagram  in 
Figure  24.  The  input  analog  signal  is  often  first  routed  through  an  analog  low  pass  anti-aliasing 
filter.  This  filtering  process  limits  the  power  spectrum  of  the  input  signal  to  half  of  the  sampling 
frequency  (Nyquist  frequency).  The  low  pass  filter  is  carefully  applied  so  as  not  to  mirror  aliased 
spectral  components  into  the  frequency  bandwidth  that  is  to  be  examined.  If  this  anti-aliasing  filter 
is  not  applied  correctly,  there  is  no  way  in  which  to  distinguish  which  part  of  the  spectrum  comes 
from  the  original  input  signal  and  which  part  of  the  spectrum  is  a  result  of  the  mirrored  spectral 
components.  An  anti-aliasing  filter  with  a  sharp-cutoff  frequency  is  required  and  can  often  be  the 
most  expensive  part  of  a  DAQ  system.  Therefore,  on  many  less  expensive  DAQ  cards,  the  sharp- 
cutoff  filter  is  simply  left  off  or  converted  to  a  more  simplified  anti-aliasing  filter  with  a  gradual 
cutoff  and  used  in  conjunction  with  oversampling  to  ensure  that  aliasing  is  not  a  problem. 
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Figure  24.  DAQ  process  block  diagrams  for  analog-digital  conversion  (top)  and  digital-analog 
conversion  (bottom). 

The  discrete  time  signal  may  then  be  subject  to  dithering,  which  adds  Gaussian  noise  to  the  signal 
with  magnitude  on  the  order  0.5  LSBrms.  This  is  used  beneficially  for  applications  involving 
averaging,  such  as  calibration  and  spectral  analysis.  Oversampling  cannot  remove  the  quantization 
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error  that  is  inherent  in  ADC,  but  by  combining  oversampling,  dithering,  and  averaging,  the  effects 
of  quantization  error  can  be  effectively  removed  and  the  stepped  nature  of  the  signal  that  results 
from  using  a  sample  and  hold  circuit  is  eliminated.  In  high-speed  DAQ  processes  that  do  not 
involve  averaging  or  spectral  analysis,  it  is  usually  a  good  idea  to  disable  dithering  to  reduce 
unwanted  noise. 

Sampling  of  the  analog  input  signal  is  achieved  using  a  “sample  and  hold”  circuit.  This 
circuit  takes  the  signal  value  at  a  particular  sampling  instant  and  holds  it  over  the  entire  sampling 
period.  The  result  is  a  stepped  signal  as  is  shown  in  Figure  25.  This  is  a  necessary  step  as  the 
analog-to-digital  converter  (ADC)  requires  a  constant  value  for  a  portion  of  the  sampling  period  so 
that  it  has  time  to  decide  on  the  correct  digital  value  to  represent  the  analog  voltage.  The  most 
common  ADC  in  high-speed  DAQ  applications  transforms  the  analog  stepped  signal  from  the 
output  of  the  sample  and  hold  circuit  into  discretely  sampled  digital  binary  values  using  a  technique 
known  as  successive  approximation.  This  process  involves  stepping  through  each  value  in  the 
binary  word  from  the  most  significant  bit  (MSB)  to  the  least  significant  bit  (LSB)  to  find  the  binary 
word  value  that  most  closely  approximates  the  numerical  analog  value.  Because  only  a  certain 
number  of  digital  values  are  available,  depending  on  the  length  of  the  binary  word  used  to  describe 
a  particular  value  (e.g.  12-bit  or  16-bit),  there  is  necessarily  a  quantization  error  during  ADC.  The 
effect  of  quantization  error  on  dimension  has  been  studied  before,  although  not  using  this  method, 
such  as  in  Kumaraswamy  et  al.  in  which  a  mathematical  relation  between  the  performance  of  a 
vector  quantizer  and  the  intrinsic  "fractal"  dimension  of  a  data  set  is  established  (Kumaraswamy  et 
al.  2004). 
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Figure  25.  Example  (idealized)  output  of  a  sample  and  hold  circuit. 
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In  many  DAQ  applications  a  digital  signal  must  be  converted  into  an  analog  signal,  e.g.  to  drive  an 
electromechanical  shaker  for  vibration  analysis  of  structures.  This,  for  example,  is  a  key  component 
of  modulated  chaotic  waves  for  use  in  ultrasonic  active  sensing  for  SHM.  The  first  step  is  a  digital- 
to-analog  conversion  (DAC),  which  generates  a  sequence  of  electrical  pulses  at  each  sampling  time 
step  whose  amplitude  is  proportional  to  the  value  of  the  corresponding  digital  samples.  The  most 
common  manner  in  which  the  amplitude  of  these  pulses  is  determined  is  to  sum  the  currents  of 
switched  current  sources  with  respect  to  an  analog  reference  voltage  and  is  known  as  a  multiplying 
DAC.  If  the  conversion  process  involves  digital  samples  with  N-bit  resolution,  the  converter 
provides  N  current  sources.  For  instance,  if  the  LSB  of  the  digital  word  is  represented  by  a  current  I, 
then  each  succeeding  source  has  current  values  of  2*1,  4*1,  all  the  way  up  to  a  value  of  2^’^*I  for  the 
MSB.  Current  sources  are  then  switched  on  or  off  depending  on  the  value  of  the  corresponding  bit 
of  the  digital  word.  Current  sources  are  preferred  over  voltage  sources  for  this  method  because  they 
have  a  higher  achievable  conversion  speed. 

Unfortunately,  in  practical  implementations  when  a  DAC  switches  from  one  voltage  level  to 
another  it  produces  glitches  in  the  output  analog  signal  due  to  released  charges.  The  largest  glitches 
occur  when  the  MSB  switches  values.  The  result  is  a  period  in  which  the  analog  output  value  is 
uncertain.  To  combat  this  problem,  a  circuit  known  as  a  "de-glitcher"  is  used  to  sample  the  output  of 
the  DAC  when  it  has  reached  a  stable  voltage  and  to  provide  an  impulse  at  this  voltage  with  a  well- 
defined  length,  usually  the  sampling  period.  This  process  is  known  as  a  zero-order  hold  and  results 
in  a  stepped  signal.  It  is  somewhat  analogous  to  the  sample  and  hold  circuit  used  in  ADC. 

If  the  stepped  analog  signal  received  from  the  output  of  the  de-glitcher  circuit  is  input  into 
an  ideal  low-pass  reconstruction  filter,  the  high  frequency  content  (steps)  of  the  signal  is  removed. 
The  low-pass  filter  interpolates  the  digital  samples  to  reproduce  the  intended  analog  output  signal. 
In  reality,  an  ideal  low-pass  filter  cannot  be  achieved,  and  the  amount  of  effort  spent  trying  to 
closely  replicate  the  ideal  filter  performance  depends  upon  the  application.  In  applications  where  the 
best  possible  reconstruction  of  the  signal  is  desired,  oversampling  techniques  which  allow  linear- 
phase  filtering  in  the  digital  domain  (FIR  filters)  are  often  used.  Again,  the  use  of  an  FIR  filter  will 
not  increase  the  dimension  of  the  signal. 
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4. 1 .2.  Data  Acquisition  Simulation 

The  simulated  DAQ  system  used  in  this  work  is  designed  to  emulate  the  DAQ  card  that  is 
being  used  in  subsequent  physical  experiments,  the  National  Instruments  PCI-6036E.  It  does  not 
incorporate  an  anti-aliasing  low  pass  filter.  Therefore,  in  this  study,  the  effects  of  a  low  pass  filter 
will  not  be  examined.  The  card  uses  standard  sample  and  hold  circuitry  and  performs  ADC  using 
successive  approximation  and  1 6-bit  resolution.  Dithering  is  automatically  enabled  and  cannot  be 
turned  off.  DAC  is  performed  using  the  multiplying  method  that  is  described  above.  No  analog  low 
pass  filters  are  available  so  these  effects  will  also  not  be  considered. 

A  Lorenz  chaotic  signal  is  created  using  the  system  of  equations  shown  in  Eq.  (42).  This  signal  has 
significant  frequency  information  from  0-10  Hz.  In  this  study,  the  simulated  "analog"  version  of  the 
signal  is  discretely  sampled  at  50  kHz.  The  effect  of  dithering  is  then  added  to  the  signal  by  adding 
noise  that  ranges  from  ±0.5*LSB.  The  sample  and  hold  circuit  is  implemented  by  sampling  the 
simulated  analog  signal  at  500  Hz.  The  ADC  is  simulated  using  a  simple  b-bit  digitize  function  that 

transforms  all  the  actual  values  of  the  signal  to  integers  between  “2  +1  and  2  .  Signals  are 

digitized  using  4,  8,  12,  16,  20,  and  24-bit  ADC.  After  each  step  in  the  DAQ  processes  the 
dimension  of  the  signal  is  verified  using  a  false  nearest  neighbors  (FNN)  approach  (Kennel  and 
Abarbanel  1992). 
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Figure  26.  Embedding  dimension  after  dithering  for  various  levels  of  ADC  resolution  (left)  and 


after  sample-and-hold  for  various  levels  of  dither  (right). 


Figure  26  (left)  shows  results  of  the  FNN  algorithm  for  the  original  "analog"  signal  as  well 
as  the  "analog"  signal  plus  dither  for  varying  levels  of  ADC  resolution.  The  signals  corresponding  to 
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4-bit  and  8-bit  ADC  resolution  are  clearly  higher  dimensional  than  the  remaining  signals.  This  is 
because  the  amplitude  of  the  dither  added  to  the  signal  is  large  and,  because  white  noise  is  infinite 
dimensional,  the  decreased  signal-to-noise  ratio  (SNR)  results  in  a  higher  dimensional  signal. 
Therefore,  using  dither  in  ADC  systems  with  low  resolution  without  averaging  is  not  suggested  and 
will  not  be  used  for  the  remainder  of  this  study.  However,  for  ADC  resolutions  of  12-bits  and 
above,  the  FNN  approach  cannot  distinguish  between  the  dithered  and  original  signals.  The  proper 
embedding  dimension  for  the  original  Lorenz  signal  is  known  to  be  /w=3  which  is  the  same  result  as 
is  seen  in  Figure  26  (left). 

The  sample  and  hold  circuit  introduces  a  stepped  nature  to  the  signal  which  can  be 
interpreted  as  high  frequency  content  at  each  of  the  discontinuities.  This  high  frequency  content  has 
an  effect  similar  to  noise  contamination  and  will  cause  the  percentage  of  FNN  to  level  off  or  even 
begin  to  increase  for  higher  embedding  dimensions.  This  effect  can  be  seen  in  Figure  26  (right).  The 
original  Lorenz  signal  and  the  dithered  signals  are  indistinguishable  after  the  sample  and  hold  is 
performed,  as  one  should  expect  since  this  was  also  the  case  after  dither  was  added.  For  this  case  it 
is  concluded  that  the  proper  embedding  dimension  is  w=3  for  all  the  signals,  despite  the  rise  in 
percent  FNN  for  higher  embedding  dimensions.  This  decision  is  made  knowing  that  there  is  some 
high  frequency  content  that  cannot  be  characterized  for  any  embedding  dimension,  and  so  the 
lowest  dimension  that  satisfactorily  unfolds  the  original  signal  content  is  used. 
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-^16-bit  Digitized  Lorenz  Signal  (dither) 
-•-20-brt  Ogitized  Lorenz  Signal  (dither) 
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Figure  27.  Embedding  dimension  after  ADC  for  different  ADC  resolutions. 
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Figure  27  shows  the  signals  after  they  have  been  put  through  the  ADC.  Clearly,  using  an 
ADC  with  low  resolution  (4‘bit,  8-bit)  results  in  increased  embedding  dimension.  This  is  a  direct 
result  of  the  large  quantization  error  for  these  low-resolution  ADC  systems.  The  discrete  time  signal 
of  quantization  errors  is  often  uncorrelated  in  complex  signals  and  can  therefore  be  characterized  as 
a  white  noise  process.  This  addition  of  white  noise  to  the  signal  reduces  the  SNR  in  much  the  same 
way  as  adding  dither  and  can  therefore  increase  the  dimension  of  the  signal.  The  SNR  decreases  by 
approximately  6  dB  for  each  bit  subtracted  from  the  digital  word  length  of  the  quantized  samples. 
Thus,  the  SNR  for  a  12-bit  quantizer  is  24  dB  higher  than  for  an  8-bit  quantizer.  For  ADC  with  a 
resolution  of  12-bits  and  above,  the  embedding  dimension  remains  unchanged. 

As  can  be  seen  below  in  Figure  28,  the  embedding  dimension  increases  as  the  frequency 
ratio  R  increases  and  the  modulation  depth  d  increases  for  a  chaotic  ultrasonic  signal  centered  at  200 
kHz.  This  result  is  expected  because  signal  bandwidth  rises  as  each  of  the  two  parameters  increase, 
as  is  seen  in  Figure  23.  As  the  bandwidth  of  the  signal  grows  it  differs  more  from  a  two- 
dimensional  sine  wave  and  becomes  increasingly  multi-dimensional.  For  chaotically  amplitude- 
modulated  signals  the  embedding  dimension  is  unaffected  during  the  DAQ  process  as  long  as  the 
card  has  12-bit  precision  or  higher. 


Figure  28.  Embedding  dimension  of  modulated  signals  using  (left)  c/=  1.0  and  (right)  R  =  0.25. 


In  order  to  verify  simulation  results  that  showed  that  there  is  no  dimension  increase  using  a  DAQ 
card  as  the  signal  goes  through  a  DAC  or  ADC  process,  a  chaotic  ultrasonic  signal  centered  at  200 
kHz  is  created  and  applied  to  the  DAC  on  a  16-bit  NI  PCI-6036E  DAQ  card.  That  signal  is  then 
converted  back  to  digital  data  using  the  DAQ  card’s  ADC.  The  estimated  embedding  dimension  of 
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the  signal  that  had  been  routed  through  the  DAC  and  ADC  was  identical  to  the  original  signal  as 
computed  by  FNN. 


4.1.3.  EfYect  of  Signal-to-Noise  Ratio 

Signal-to-noise  ratio  (SNR)  can  be  defined  in  decibel  form  using  the  following  Eq.  (43): 


(43) 


where  A  is  the  root-mean-square  (RMS)  value  of  the  signal  or  noise  amplitude.  Therefore  an  order 
of  magnitude  increase  in  signal  amplitude  will  result  in  the  SNR  rising  by  20  dB.  To  identify  the 


effect  of  SNR  on  signal  dimension  white  noise  is  added  to  the  standard  Lorenz  signal  and  then  to  a 


chaotic  ultrasonic  signal  centered  at  200  kHz.  As  can  be  seen  in  Figure  29,  embedding  dimension 
decreases  as  SNR  increases.  This  agrees  with  prior  observations  made  about  the  effect  of  SNR  on 
embedding  dimension. 


Figure  29.  Embedding  dimension  of  (left)  standard  Lorenz  signal  and  (right)  a  200  kHz 
chaotic  amplitude-modulated  signal  for  various  levels  of  SNR. 

This  section  outlined  a  method  for  creating  chaotic  amplitude-modulated  ultrasonic  signals. 
Several  signal  creation  parameters  are  introduced  and  their  effects  on  signal  characteristics  such  as 
power  spectral  density  are  detailed.  A  study  of  the  embedding  dimension  of  chaotic  ultrasonic 
signals  (through  the  creation  and  data  acquisition  processes)  is  also  presented.  The  concepts  of  state- 
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space  embedding  and  attractor  reconstruction  are  introduced  as  a  means  to  estimate  embedding 
dimension.  A  low  signal  dimension  is  necessary  to  provide  reduced  processing  and  storage 
requirements  and  increase  the  ability  of  the  system  to  resolve  damage  via  statistically  robust 
computation  of  features.  It  is  shown  that  the  DAQ  process  does  not  substantially  alter  the  embedded 
dimension  of  a  chaotic  signal  as  long  as  the  ADC/D  AC  resolution  is  at  least  12-bits. 

4.2.  Damage  Detection/Classification  Scheme 

Fundamental  time  series  modeling  where  prediction  power  is  required  is  often  tackled  via 
auto-regressive  (AR)  models,  which  have  previously  been  shown  to  be  useful  as  a  damage  classifier 
(Sohn  and  Farrar  2001;  Sohn  et  al.  2001).  The  discretely  observed  output  time  series  x(a7)  is 
modeled  with  an  AR  model  (e.g.,  a  finite  impulse  response  filter)  of  the  form: 

p 

xf'j)  =  -  /)  +  e(A7)  (44) 


where  p  is  the  order  of  the  AR  model  with  associated  coefficients  and  residual  error  e(A7).  This 
dissertation  examines  various  model  orders  p  to  determine  the  most  suitable  model  order  for  each 
particular  damage  detection  scheme.  The  AR  coefficients  are  estimated  through  minimization  of  the 
sum-of-squared  forward  prediction  errors  (Brockwell  and  Davis  1991).  All  signals  are  normalized 
through  division  of  the  standard  deviation  of  the  signal  before  use  of  the  AR  model. 

The  calculated  AR  coefficients  themselves  can  be  used  a  damage-sensitive  feature.  This 
dissertation  examines  the  ability  of  a  two-part  classification  scheme  to  determine  the  existence  and 
size  of  damage  in  a  structure  using  a  supervised  learning  paradigm.  This  classification  scheme  can 
also  be  employed  with  an  unsupervised  learning  technique  to  classify  whether  a  structural 
connection  for  which  there  is  no  baseline  data  is  damaged  or  undamaged.  The  classification 
approach  uses  features  derived  from  an  autoregressive  (AR)  model  coefficient  vector  cross¬ 
assurance  criterion  as  well  a  group  classification  scheme  with  its  basis  in  information  theory. 
Throughout  this  section,  data  from  a  physical  experiment  on  several  composite  adhesively-bonded 
wing-to-spar  test  structures,  that  will  be  discussed  in  detail  later  in  this  report,  is  used  to  illustrate  the 
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effectiveness  and  necessity  of  the  two-part  damage  classification  scheme.  Two  of  these  structures 
have  various  disbonds  built  into  the  bond  line  through  the  use  of  Teflon  inserts  and  one  of  the 
structures  has  no  damage  built  in.  The  experiment  outlined  later  uses  the  two-part  supervised 
learning  classification  scheme  outlined  in  this  section  to  not  only  to  classify  disbond  size  but  also  to 
classify  whether  a  bond  for  which  there  is  no  baseline  data  is  undamaged  or  has  some  form  of 
disbond. 

4.2. 1 .  Vector  consistency  criterion 

The  first  step  in  the  damage  classification  algorithm  is  based  on  a  comparison  of  AR 
coefficient  vectors.  Figure  30  shows  averaged  AR  coefficients  that  are  acquired  using  an  AR  model 
order  of  p=\0  for  three  different  bond  conditions,  including  an  undamaged  bond  as  well  as  two 
different  disbond  sizes  on  a  single  adhesively-bonded  composite  test  structure.  Figure  30  (right) 
shows  averaged  AR  coefficients  that  are  acquired  using  an  AR  model  order  of  p=l0  for  three 
different  undamaged  bond  conditions  on  a  similarly  configured  test  structure.  The  figures  show  a 
clear  difference  in  the  shape  of  the  AR  coefficient  vectors  for  different  bond  conditions  whereas  all 
undamaged  bond  conditions  result  in  similarly  shaped  AR  coefficient  vectors.  This  result  suggests 
that  a  metric  which  compares  the  shape  of  two  AR  coefficient  vectors  will  be  a  damage-sensitive 
feature. 


Figure  30.  Averaged  AR  coefficient  vector  comparison  for  (left)  different  bond  conditions  and 
(right)  three  undamaged  bond  conditions. 
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A  feature  based  on  the  modal  assurance  criterion  (MAC)  that  is  commonly  used  in  modal 
analysis  to  provide  a  measure  of  consistency  between  estimates  of  a  modal  vector  was  proposed 
(Allemang  2003) .  The  MAC  is  a  scalar  value  detailing  the  linear  relationship  between  a  test  modal 
vector  and  a  reference  modal  vector  as  follows: 


MAC„ 


W.YW.}{v.rM 


(45) 


where  \i/t  is  the  test  modal  vector,  xf/r  is  the  reference  modal  vector,  and  H  indicates  the  Hermitian 
operator  or  complex  conjugate.  The  modal  assurance  criterion  can  have  values  ranging  from  zero, 
representing  no  consistent  correspondence,  to  unity,  representing  complete  consistent 
correspondence.  Therefore  if  the  two  vectors  are  consistently  related,  then  the  MAC  value  should 
approach  unity.  This  study  uses  a  vector  consistency  criterion  (VCC)  similar  to  the  MAC  in  order 
to  measure  the  orthogonality  of  AR  coefficient  vectors.  In  this  case  the  Hermitian  operator 
simplifies  to  a  vector  transpose  because  the  values  of  the  AR  coefficient  vector  are  always  real.  The 
fact  that  the  VCC  feature  can  only  have  a  value  between  zero  and  unity  makes  it  useful  as  an 
extracted  feature  because  scaling  effects  of  the  original  signal  will  not  affect  any  decision  that  needs 
to  be  made  regarding  the  bond  condition  assessment.  The  hypothesis  being  tested  is  that  two  AR 
coefficient  vectors  should  be  very  consistent  (nearly  orthogonal)  if  they  are  derived  from  modeled 
responses  of  structural  connections  with  the  same  level  (or  lack  of)  damage. 


4.2.2.  Step  one:  damage  existence  classification 

The  determination  of  whether  a  connection  of  unknown  condition  is  undamaged  or  damaged 
is  accomplished  using  the  following  method.  First,  a  set  of  distinct  250-microsecond-long  input 
signals  are  created  from  the  data-generating  process  that  has  been  previously  described.  For  each  of 
these  input  signals  a  structural  response  is  recorded  once  under  each  of  the  known  undamaged  and 
damaged  conditions  on  a  particular  specimen,  using  a  different  actuator/sensor  pair  for  each  bond, 
as  well  as  a  second  time  for  the  undamaged  condition.  AR  coefficients  are  then  estimated  for  each 
of  these  structural  responses.  These  sets  of  AR  coefficients  form  a  database  of  structural  responses 
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for  known  structural  conditions  with  each  input  signal  that  will  be  used  to  determine  the  condition 
of  an  unknown  connection  in  a  supervised  learning  fashion. 

For  the  first  step  of  the  bond  condition  assessment  algorithm,  which  will  identify  only 
whether  an  unknown  connection  is  damaged  or  undamaged,  a  distribution  of  VCC  values  is  then 
calculated.  This  is  accomplished  by  comparing,  using  Eq.  (45),  each  of  the  AR  coefficient  vectors 
from  the  first  set  of  waveforms  to  probe  the  undamaged  connection  with  each  of  the  AR  coefficient 
vectors  from  the  second  set  of  waveforms  to  probe  the  undamaged  connection.  This  calculation 
results  in  a  large  distribution  of  VCC  values  that  describe  the  auto-correlation  between  the  two  sets 
of  AR  coefficient  vectors  that  are  both  derived  from  the  undamaged  connection  condition.  A  new 
set  of  input  signals  (created  from  the  same  underlying  process  as  the  supervised  learning  database 
input  signals)  are  then  applied  to  the  structure  for  a  structural  connection  whose  damage  state  needs 
to  be  determined.  AR  coefficient  vectors  are  again  calculated  for  each  of  these  new  structural 
responses.  The  AR  coefficient  vectors  from  the  original  undamaged  connection  condition  are  then 
compared  with  the  new  set  of  AR  coefficient  vectors  from  the  unknown  connection  condition  by 
calculating  VCC  values.  This  computation  results  in  a  second  distribution  of  VCC  values  that 
describe  the  correlation  between  the  AR  coefficient  vectors  of  the  undamaged  connection  condition 
and  the  AR  coefficient  vectors  of  the  unknown  connection  condition. 

A  statistical  measure  must  now  be  devised  that  can  properly  classify  a  damaged  structural 
connection  and  at  the  same  time  properly  classify  an  undamaged  connection  in  a  different  geometric 
location  using  a  different  MFC  sensor/actuator  pair  than  the  baseline  undamaged  connection.  Figure 
31  (left)  shows  empirical  probability  density  functions  (pdfs)  of  VCC  values  for  three  different  bond 
conditions,  including  an  undamaged  bond  as  well  as  two  different  disbond  sizes  from  a  single 
adhesively-bonded  composite  test  structure.  Figure  31  (right)  shows  empirical  pdfs  of  VCC  values 
for  three  different  undamaged  bond  conditions  in  different  geometric  locations  on  a  different,  but 
similarly  configured,  test  structure.  It  is  clear  that  while  there  are  differences  in  the  pdfs  of  the 
undamaged  bonds  in  Figure  31  (left),  the  difference  between  the  pdfs  of  the  undamaged  bond  and 
two  disbond  sizes  in  Figure  31  (right)  are  much  greater.  A  statistical  decision  boundary  that  can 
account  for  the  small  differences  in  the  pdfs  due  to  geometric  conditions  as  well  as  MFC  bonding 
characteristics  while  still  correctly  classifying  the  damaged  bonds  would  be  beneficial.  The  pdfs  are 
not  normally  distributed  as  can  be  seen  from  the  elongated  left  tail.  Therefore  standard  statistical 
tests  that  assume  normality  such  as  the  Student’s  t-test  cannot  be  used  in  this  case.  There  are  several 
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non-parametric  tests  that  can  be  employed  for  non-Gaussian  distributions  such  as  the  Kolmogorov- 
Smirnov  test  or  the  Mann- Whitney  U  test,  but  for  this  application  these  tests  decidedly  classify  the 
various  undamaged  pdfs  as  being  from  different  parent  distributions.  A  solution  to  this  problem  is  to 
set  a  lower  confidence  limit  on  an  undamaged  pdf  so  that  the  other  undamaged  pdfs  would  not 
result  in  a  significant  number  of  outliers  (VCC  values  less  than  the  lower  bound)  but  still  produce  a 
significant  number  of  outliers  for  the  damaged  pdfs. 


VCC  Value  VCC  Value 


Figure  31.  VCC  value  empirical  pdf  comparison  for  (left)  different  bond  conditions  and  (right) 
three  undamaged  bond  conditions. 

Extreme  value  statistics  (EVS)  are  used  to  set  the  lower  confidence  limit  on  the  non -normal 
distributions  that  are  produced  using  the  above  outlined  method.  EVS  is  used  in  this  analysis  to 
accurately  model  the  behavior  of  the  feature  distribution’s  tails.  The  basis  of  this  branch  of  statistics 
stems  from  the  following  situation.  If  a  moving  window  is  taken  along  a  vector  of  samples  and  an 
extremum  value  (minimum  or  maximum)  is  selected  from  each  of  these  windows,  the  induced 
cumulative  density  function  of  the  extrema  of  the  samples,  as  the  number  of  vector  samples  tends  to 
infinity,  asymptotically  converges  to  one  of  three  possible  distributions:  Gumbel,  Weibull  or 
Frechet  (Castillo  1988).  The  following  equations  detail  the  Gumbel,  Weibull  and  Frechet 
distributions  for  minima  because  that  is  what  this  study  uses. 
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where  ^ ^ ,  and  ^  are  the  model  parameters  that  are  estimated  from  the  data.  There  are  similar 
functions  that  apply  to  the  converged  distributions  for  the  maxima  of  the  sample  sets. 

The  appropriate  distribution  is  chosen  by  plotting  the  cumulative  distribution  function  (cdf) 
of  the  extracted  vector  of  minima  on  the  probability  paper  for  a  Gumbel  distribution  (Castillo  1988). 
Probability  paper  is  designed  so  that  the  cdf  for  a  particular  cumulative  distribution  will  plot  in  a 
linear  fashion  by  transforming  the  vertical  scale  of  the  plot.  In  this  case,  if  the  vector  of  minima  has 
a  Gumbel  minimum  distribution,  the  cdf  will  be  a  straight  line  when  plotted  on  Gumbel  probability 
paper.  Otherwise,  the  cdf  will  have  an  associated  curvature.  If  this  curvature  is  concave,  the  feature 
vector  has  a  Weibull  minimum  distribution.  Similarly,  if  the  curvature  is  convex  the  feature  vector 
has  a  Frechet  minimum  distribution.  In  this  study  the  distribution  of  minima  taken  from  the  parent 
undamaged  VCC  distributions  has  a  Gumbel  minimum  distribution.  Model  parameters  are  then 
estimated  by  fitting  the  chosen  distribution  to  the  data  as  outlined  in  Castillo  (Castillo  1988). 

Once  the  model  parameters  are  chosen,  it  is  possible  to  generate  confidence  limits  that  can 
be  applied  to  the  distribution.  These  limits  are  far  more  accurate  than  those  obtained  when  assuming 
a  Gaussian  distribution.  The  threshold  corresponding  to  a  specific  confidence  level  for  the  Gumbel 
minimum  distribution  is  given  by  the  following  equation  (Worden  et  al.  2002): 


Gumbel:  x_:_  =  A  +  5  In 
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where  n  is  the  window  size  used  to  extract  the  maxima,  and  ^  is  the  associated  Type  I  error  of  the 
confidence  limit.  Type  I  error  is  simply  the  percentage  of  false  positives  that  are  expected  to  appear. 
For  example,  when  using  a  95%  confidence  interval  the  Type  I  error  is  expected  to  be  5%  ( 
a  =  0.05 )  3y  using  the  confidence  limit  calculated  through  EVS  the  damage  state  classification 
algorithm  is  able  to  differentiate  between  undamaged  and  damaged  connections.  To  further 
differentiate  between  different  sizes,  levels  or  types  of  damage  a  second  classification  scheme  is 
used. 

4.2.3.  Step  two:  damage  size/tvpe  classification 

This  part  of  the  project  developed  a  novel  statistical  classification  technique  with  its  basis  in 
information  theory  to  classify  different  sizes,  levels  or  types  of  damage  in  bonds.  The  classification 
method  is  based  on  a  fundamental  theorem  of  Shannon’s  information  theoiy  that  states  that  “the 
best  compression  for  any  given  data  set  comes  from  a  codebook  designed  exactly  for  the  statistics  of 
that  source;  any  other  codebook  will  give  worse  results”  (Shannon  1949).  Recently,  information 
theoiy  has  been  applied  to  continuous  time  signals,  where  compression  performance  is  related  to 
prediction  error,  and  a  codebook  is  the  model  for  a  source  that  produces  time  series  data.  In  this 
study  the  source  is  a  chaotic  guided  wave  propagating  through  a  structural  connection.  These  time 
series  can  then  be  classified  using  prediction  error  as  a  means  of  virtual  data  compression  (via  the 
AR  coefficient  vectors).  The  model  for  the  response  time  series  is  the  same  auto-regressive  model 
used  for  the  undamaged/damaged  discrimination  that  has  already  taken  place.  This  idea  leads  to  a 
procedure  for  classifying  time  series  using  cross-prediction  error,  as  literal  “data  compression”  in 
the  sense  of  maximal  information  recovery  on  decompression  is  not  actually  necessaiy,  just  its 
“virtual”  performance  as  a  representative  data  model.  Of  course,  the  better  the  underlying  statistical 
model  is  the  more  the  classification  performance  will  improve. 

The  entire  damage  state  classification  technique  can  be  summarized  as  follows.  The  same 
database  of  AR  coefficients  that  was  created  using  a  set  of  chaotic  ultrasonic  guided  waves  is  again 
used  to  decide  the  state  of  the  damaged  connection.  In  this  case  only  AR  coefficients  from  the 
damaged  paths  are  used  because  it  is  now  known  that  the  test  path  has  some  form  of  damage  after 
using  the  first  classification  technique.  A  new  input  signal  (created  from  the  same  underlying 
process  as  the  training  database  input  signals)  is  then  applied  to  the  structure  along  the  known 
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damaged  bond  condition  path.  One  of  the  sets  of  AR  coefficients  in  the  training  database  for  each  of 
the  paths  with  disbond  damage  is  then  used  to  estimate  the  structural  response  to  the  new  input 
signal.  One  set  of  coefficients  from  the  training  database  will  minimize  the  sum  of  the  squared 
residual  errors  and  these  coefficients  are  associated  with  a  known  size,  level  or  type  of  damage  in 
the  training  database.  This  damage  state  is  then  classified  as  the  state  of  the  unknown  connection 
that  is  being  probed.  This  comparison  takes  place  for  each  of  the  remaining  input  signals  in  the 
training  database.  This  entire  process  is  then  repeated  using  a  large  set  of  input  signals  that  are 
imparted  to  the  structure  in  its  unknown  damage  state  (the  same  data  as  was  used  in  the  first 
classification  step).  The  votes  for  each  damage  state  are  then  summed  and  the  plurality  of  votes  is 
the  estimated  condition  of  the  bond.  The  entire  two-part  statistical  damage  state  classification 
paradigm  is  depicted  visually  Figure  32. 


Figure  32.  Ultrasonic  statistical  classification  paradigm. 


4.3.  Experimental  Validation  on  Adhesively-Bonded  Composite  Joints 

4.3. 1 .  Single  Composite  Plate  Tests 

This  study  applies  the  method  outlined  above  to  classify  different  bond  state  damage 
conditions  of  an  adhesively-bonded  composite  joint,  including  various  disbond  sizes  and  a  poorly- 
cured  bond.  The  test  structure  is  analogous  to  a  wing  skin-to-spar  bonded  Joint  and  consists  of  a 
custom  made  8-ply  [0/±45/0]s  carbon  fiber-reinforced  plastic  (CFRP)  plate  that  is  bonded  to  a 
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tubular  CFRP  spar  (see  Figure  10).  The  CFRP  plate  measures  12  in  x  12  in  x  0.04  in,  and  the  spar 
has  an  outer  square  diameter  of  2  in  x  2  in  with  a  wall  thickness  of  0.25  in.  The  bonded  area 
contains  two  different  sizes  of  disbond  (0,25  in  and  1 ,0  in  )  created  using  Teflon  inserts  as  well  as  a 
portion  of  the  bond  that  has  been  poorly  cured  by  using  only  50%  of  the  prescribed  level  of 
hardener.  An  active  excitation  signal  is  imparted  to  the  structure  through  a  MFC  patch  on  one  side 
of  the  bonded  joint  and  sensed  using  an  equivalent  MFC  patch  on  the  opposite  side  of  the  joint  and 
there  is  an  MFC  actuator/sensor  pair  for  each  bond  condition  to  be  identified. 

Two  MFC  actuator/sensors  were  attached  to  the  structure  on  each  side  of  the  bond  for  each 
of  the  three  damaged  bond  states  as  well  as  one  pair  for  the  undamaged  bond  area  for  a  total  of  eight 
MFC  patches.  Each  MFC  patch  is  affixed  to  the  structure  one  inch  from  the  spar  bond  line.  Each 
input  signal  is  applied  to  the  structure  25  times  and  then  averaged  and  filtered  to  reduce 
experimental  noise.  The  actuation  signal  is  created  by  the  output  channel  of  a  National  Instruments 
PCl-61 10  DAQ  card  at  a  rate  of  4  MHz  and  routed  through  a  Krohn-Hite  7602  wideband  power 
amplifier.  This  amplified  signal  is  sent  to  the  actuation  MFC  while  the  sensing  MFC  simultaneously 
samples  the  structural  response  at  a  rate  of  4  MHz. 

For  the  following  experiments,  30  time  response  histories  are  created  using  a  specific  set  of 
signal  creation  parameters.  Half  of  these  30  generated  response  time  histories  are  selected  as 
database  training  inputs.  The  remaining  15  responses  are  used  as  test  set  inputs.  All  four  bond 
conditions  are  probed  three  times  with  the  first  two  time  history  responses  being  used  for  the 
training  database  and  the  remaining  responses  used  to  determine  bond  condition  for  the  test  set 
inputs.  In  each  case  the  actual  bond  condition  is  determined  based  solely  on  knowledge  it  acquired 
from  the  15  database  training  inputs. 

In  order  to  determine  the  set  of  input  signal  parameters  that  produce  optimal  bond  state 
identification,  the  statistical  classification  scheme  is  employed  for  an  array  of  signal  parameters  and 
the  overall  ability  of  each  set  of  parameters  to  correctly  identify  each  bond  state  is  then  examined. 
Previous  numerical  and  experimental  work  done  by  Matt,  Bartoli  and  Lanza  di  Scalea  on  the  same 
test  structure  showed  that,  when  using  through  transmission  characteristics  of  simple  toneburst 
inputs,  an  excitation  frequency  of  200  kHz  was  best  able  to  identify  bond  state  condition  due  to  a 
coupling  of  So  and  Ai  modes  at  that  frequency  (Matt  et  al.  2005;  Lanza  di  Scalea  et  al.  2007).  In  this 
study  a  number  of  excitation  center  frequencies  from  100  kHz  to  300  kHz  are  considered,  as  well  as 
frequency  ratios  ranging  from  0.01  (nearly  sinusoidal)  to  0.50  (highly  chaotic).  When  considering 
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the  feature  extraction  parameters  a  number  of  AR  model  orders  from  3  to  25  are  employed.  Figure 
33  depicts  the  ability  of  the  statistical  classification  scheme  to  correctly  identify  bond  condition  by 
plotting  the  percentage  of  correct  votes  against  the  various  signal  input  and  feature  extraction 
parameters  already  discussed. 
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Figure  33.  Classification  scheme  results  showing  percentage  of  total  correct  votes  vs.  (top  left) 
frequency  ratio,  (top  right)  AR  model  order,  (lower  left)  center  frequency  and  (lower  right)  center 
frequency  for  each  individual  bond  condition. 

The  top  left  plot  of  Figure  33  shows  the  ability  of  the  statistical  classification  scheme 
employed  together  with  chaotic  insonification  to  correctly  identify  bond  condition  for  various 
frequency  ratios.  As  previously  mentioned,  the  lower  the  frequency  ratio  the  more  sinusoidal  the 
input  signal  appears  and  conversely  the  higher  the  frequency  ratio  the  more  chaotic  the  input  signal 
appears.  This  signal  parameter  characteristic  occurs  because  as  the  frequency  ratio  becomes  small, 
the  chaos  operates  on  a  longer  time  scale.  If  a  small  enough  time  window  is  used  to  create  the 
signal,  a  very  small  R  value  will  result  in  the  chaotic  time  scale  being  much  longer  than  that  of  the 
entire  signal.  In  this  case  the  signal  will  appear  to  be  almost  sinusoidal.  The  opposite  effect  holds 
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true  if  a  large  frequency  ratio  is  used  because  the  time  scale  the  chaos  operates  on  will  be  much 
smaller  than  the  length  of  the  signal  and  the  signal  will  therefore  contain  more  chaotic  information. 
It  is  clear  that  a  signal  with  more  chaotic  information,  when  paired  with  this  pattern  recognition 
scheme,  is  best  able  to  determine  true  bond  state  condition.  There  is  negligible  difference  between  a 
frequency  ratio  of  0.3  and  0.5  but  the  overall  trend  that  a  more  chaotic  signal  provides  greater 
classification  ability  is  a  sufficient  conclusion.  The  upper  right  plot  depicts  the  effect  of  AR  model 
order  on  the  damage  identification  ability  of  the  statistical  classification  scheme.  In  general,  if  too 
low  of  an  AR  model  is  used  (anything  less  than  10)  the  percentage  of  total  correct  votes  has  a 
noticeable  decrease.  When  averaged  over  all  frequency  ratios  and  center  frequencies  an  AR  model 
order  of  10  does  slightly  better  than  higher  AR  model  orders.  However,  the  effectiveness  of  a 
particular  model  order  is  clearly  dependent  on  the  particular  frequencies  in  which  the  chaotic 
information  is  present.  A  closer  investigation  of  the  data  reveal  that,  depending  on  the  signal  input 
parameters,  an  AR  model  order  somewhere  between  10  and  20  will  produce  optimal  classification 
of  each  particular  bond  state. 

The  lower  left  plot  in  Figure  33  illustrates  the  effect  of  various  center  excitation  frequencies 
on  the  ability  to  correctly  classify  bond  condition.  The  graph  shows  that  center  frequency  is  a  highly 
important  parameter  as  there  appears  to  be  “sweet  spots”  where  the  algorithm  correctly  classifies 
almost  100%  of  individual  test  cases  (150  IcHz  and  200  kHz)  as  well  as  several  frequencies  which 
do  not  perform  as  well  including  175  kHz  which  is  intermediate  to  the  two  sweet  spots  in  the 
frequency  regime.  As  previously  mentioned,  numerical  simulations  (Matt  et  al.  2005)  showed  a 
mode  coupling  occurring  at  200  kHz  which  could  be  the  reason  for  the  excellent  performance  of  the 
damage  detection  algorithm  across  all  frequency  ratios  and  model  orders  at  that  frequency.  The 
same  study  does  not  indicate  a  mode  coupling  effect  at  1 50  kHz.  However,  these  are  only  numerical 
results,  and  the  true  physical  structure  may  contain  modal  coupling  effects  at  150  kHz.  Finally,  the 
lower  right  plot  of  Figure  33  shows  that  the  undamaged  bond  condition  is  correctly  identified  for 
almost  100%  of  the  individual  input  signals  for  any  center  frequency.  The  ability  to  discern  the 
poorly  cured  bond  state  appears  to  decrease  as  center  frequency  increases  and  the  small  disbond 
case  is  the  most  often  misidentified  bond  condition,  as  would  be  expected  due  to  its  size. 

4.3.2.  Temperature  Variability 
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The  experimental  apparatus  used  in  this  study  is  meant  to  simulate  a  wing  skin-to-spar 
bonded  joint  that  is  used  in  aerospace  applications.  Therefore,  if  this  statistical  classification  scheme 
using  chaotic  insonification  is  to  be  used  for  in  situ  health  monitoring,  including  in-flight  and  on  the 
ground,  its  behavior  and  effectiveness  relative  to  external  temperature  must  be  considered.  This 
study  therefore  employs  the  use  of  a  thermal  chamber  to  test  the  efficacy  of  the  above  outlined 
method  for  temperatures  varying  from  -40°C  to  40°C,  a  range  that  covers  most  flight  and  ground 
conditions.  For  the  thermal  chamber  test  a  center  frequency  of  200  kHz,  frequency  ratio  of  0.5  and 
AR  model  order  of  10  are  employed  as  these  parameters  produced  near  100%  correct  identification 
of  individual  test  signals  for  the  room  temperature  experiments.  Figure  34  shows  the  initial  thermal 
chamber  results  for  cases  in  which  the  baseline  AR  coefficients  are  determined  at  various 
temperatures.  It  is  clear  that  the  damage  detection  algorithm  only  performs  well  at  temperatures 
close  to  that  of  the  temperature  at  which  the  baseline  AR  coefficients  are  determined.  A  temperature 
change  of  more  than  10  degrees  Celsius  usually  results  in  significant  reduction  of  classification 
effectiveness.  However,  when  the  proper  baseline  AR  coefficients  are  used  for  a  particular 
temperature  the  statistical  classification  method  still  achieves  near  100%  correct  classification  of 
individual  test  signals. 


Figure  34.  Percentage  of  total  correct  votes  vs.  temperature  for  cases  in  which  the  baseline  AR 
coefficients  are  determined  at  various  temperatures. 

It  is  therefore  possible  to  correct  for  temperature  variations  by  storing  baseline  AR  coefficients  for 
the  range  of  temperatures  that  the  intended  apparatus  is  likely  to  see  during  real  world  applications. 
It  still  must  be  determined  how  many  sets  of  baseline  AR  coefficients  must  be  stored  and  what  the 
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proper  method  for  interpolation  between  these  temperatures  is.  Therefore,  this  study  examines  using 
baseline  AR  coefficients  taken  every  20  degrees  Celsius  (■40°C,  -20°C,  0°C,  20°C,  40°C)  which  is 
referred  to  in  the  following  figures  as  Stempbase  and  using  baseline  AR  coefficients  taken  every  10 
degrees  Celsius  (-40®C,  -30®C,  -20®C,  -lO^C,  0®C,  10®C,  20®C,  30®C,  40°C)  which  is  referred  to  in 
the  following  figures  as  9tempbase.  Testing  was  carried  out  at  5  degree  Celsius  temperature 
increments.  It  can  be  seen  in  Figure  35  that  the  value  of  each  AR  coefficient  tends  to  shift  in  a  fairly 
reliable  manner,  usually  in  a  near  linear  fashion,  as  the  temperature  changes.  It  is  then  possible  to 
create  a  baseline  of  AR  coefficients  at  all  temperatures  using  polynomial  fits  from  testing  carried 
out  at  only  a  few  intermediate  temperatures. 
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Figure  35.  Actual  coefficient  values  vs.  temperature  for  the  first  (left)  and  second  (right)  AR 
coefficients. 

Figure  36  shows  the  effectiveness  of  using  various  order  polynomial  fits  (from  linear  to 
fourth  order)  to  determine  baseline  AR  coefficients  for  intermediate  temperatures.  It  is  clear  that 
using  baseline  AR  coefficients  determined  every  10  degrees  Celsius  produces  better  results  than 
using  baseline  AR  coefficients  determined  every  20  degrees  Celsius.  In  both  cases  a  fourth  order 
polynomial  produces  the  highest  percentage  of  total  correct  votes  over  all  temperatures  with  no 
lower  result  than  92%  for  the  Stempbase  case  and  no  lower  result  than  98%  for  the  9tempbase  case. 
However,  because  the  AR  coefficients  were  being  fit  to  a  polynomial  the  optimal  AR  coefficients 
for  a  particular  temperature  were  not  always  used,  as  the  polynomial  does  not  necessarily  go 
directly  through  this  known  data  point. 
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Figure  36.  Percentage  of  total  correct  votes  for  various  fitting  techniques  using  the  Stempbase  (a, 
left)  and  9tempbase  (b,  right)  methods. 


Therefore,  a  lookup  table  method  in  which  the  AR  coefficients  for  temperatures  in  the  baseline 
database  are  directly  used,  and  any  intermediate  temperatures  are  simply  linearly  interpolated,  is 
also  examined.  Figure  37  shows  the  result  of  using  this  method.  It  is  clear  that  even  in  the 
Stempbase  case  (baseline  AR  coefficients  determined  in  20  degree  Celsius  increments)  that  linear 
interpolation  of  AR  coefficients  for  intermediate  temperatures  results  in  100%  correct  classification 
of  individual  test  cases.  This  result  is  not  surprising  given  Figure  35  shows  that  the  AR  coefficients 
trend  almost  piece-wise  linear  as  temperature  changes.  This  study  therefore  concludes  that  the  best 
method  for  temperature  correction  is  to  store  baseline  coefficients  for  a  small  number  of 
temperatures  that  span  the  known  operational  values  of  the  structure  and  then  to  linearly  interpolate 
via  a  lookup  table  the  AR  coefficients  necessary  to  correctly  classify  bond  condition  for  any 
intermediate  temperatures. 
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Figure  37.  Percentage  of  total  correct  votes  for  direct  linear  lookup  table  using  the  Stempbase  and 
9tempbase  methods. 

4.3.3.  Multiple  Composite  Plates 

The  previous  section  showed  the  ability  of  the  bond  state  assessment  algorithm  outlined  in 
this  study  to  compensate  for  temperature  variations  (variability  due  to  the  environment).  Other  types 
of  variability  that  can  affect  the  outcome  of  this  scheme  include  unintentional  inhomogeneities  in  a 
single  manufactured  specimen  (“within-unit”  variability)  and  inconsistencies  in  the  manufacturing 
process  itself  across  specimen  lots  (“unit-to-unit”  variability).  Within-unit  variability  includes  the 
effects  of  manufacturing  variability  along  the  bond  line,  MFC  bond  condition  inconsistency,  and 
geometric  effects  that  result  from  different  placements  of  the  MFC  sensor/actuator  pairs.  This  study 
examines  the  ability  of  the  classification  scheme  to  not  only  to  classify  disbond  size  but  also  to 
compensate  for  the  different  sources  of  variability  discussed  above  (such  as  specimen-to-specimen). 
This  classification  is  performed  using  data  from  a  similarly  configured  composite  bond  for  which 
baseline  data  are  available.  The  test  structures  are  again  analogous  to  a  wing  skin-to-spar  bonded 
joint  and  there  is  an  MFC  actuator/sensor  pair  for  each  bond  condition  to  be  identified  on  three 
composite  wing-to-spar  adhesive  bonded  Joints.  The  test  structures  consist  of  a  carbon  fiber- 
reinforced  plastic  (CFRP)  plate  (manufactured  by  McMaster-Carr)  that  is  bonded  to  a  tubular  CFRP 
spar  using  Loctite®  Hysol®  9462  Two  Component  Epoxy  Adhesive*.  The  CFRP  plate  measures 
12in  X  12in  x  0.016in  square  and  the  spar  again  has  an  outer  square  diameter  of  2in  x  2in  with  a  wall 
thickness  of  0.25in.  The  bonded  area  for  two  of  the  specimens  contains  two  different  sizes  of 
disbond  (0.25  in^  and  1 .0  in^)  that  are  created  using  Teflon  inserts.  The  structures  with  disbonds  are 
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labeled  Plate  1  and  Plate  2.  The  third  test  structure  (Plate  3)  was  manufactured  with  no  disbonds  so 
all  paths  are  undamaged  and  geometric  effects  can  be  studied.  The  locations  of  the  MFC 
actuators/sensors  for  Plate  1  and  Plate  2  are  pictured  in  Figure  38.  C-scans  that  depict  the  bond  state 
of  each  specimen  can  be  seen  in  Figure  39. 
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Figure  38.  Experimental  platform  showing  location  of  bottom-mounted  MFC  actuators  on  Plate  1 
(left)  and  top-mounted  MFC  actuators  on  Plate  2  (right). 


Figure  39.  C-scan  images  showing  locations  and  sizes  of  disbonds  for  Plate  1  (top)  and  Plate  2 
(middle)  as  well  as  undamaged  Plate  3  (bottom). 


Plate  1  features  two  MFC  actuator/sensors  that  are  attached  to  the  CFRP  plate  on  each  side 
of  the  bond  for  both  of  the  damaged  bond  states  as  well  as  one  pair  for  the  undamaged  bond  area  for 
a  total  of  six  MFC  patches.  The  MFC  patches  are  located  on  the  spar  side  of  the  bond  as  seen  in 
Figure  38.  The  MFC  patches  on  Plate  2  are  located  on  the  opposite  side  of  the  CFRP  plate  from  the 
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bonded  area.  Plate  3  also  has  six  MFC  patches  that  are  bonded  on  the  spar  side  of  the  structure,  as  in 
Plate  1.  Each  of  the  MFC  patches  are  bonded  to  the  structure  using  Loctite®  Hysol®  E120HP™ 
epoxy  adhesive.  Each  MFC  patch  is  affixed  to  the  structure  one  inch  from  the  spar  bond  line.  Each 
input  signal  is  applied  to  the  structure  25  times  and  then  averaged  and  filtered  to  reduce 
experimental  noise. 

For  the  following  experiments,  in  order  to  correctly  identify  bond  condition,  90  time 
response  histories  are  created  using  a  specific  set  of  signal  creation  parameters.  Half  of  these  90 
generated  response  time  histories  are  selected  as  database  training  inputs.  The  remaining  45 
responses  are  used  as  test  set  inputs.  All  three  bond  conditions  on  each  structure  are  probed  three 
times  with  the  first  two  time  history  responses  being  used  for  the  training  database  and  the 
remaining  responses  used  to  determine  bond  condition  for  the  test  set  inputs.  In  each  case  the  actual 
bond  condition  is  determined  based  solely  on  knowledge  it  acquired  from  the  45  database  training 
inputs. 

A  preliminary  study  is  undertaken  to  determine  the  set  of  input  signal  parameters  and 
extracted  feature  parameters  (such  as  AR  model  order)  that  produce  optimal  bond  state 
identification  for  the  new  set  of  test  structures.  This  study  uses  receiver  operating  cheiracteristic 
(ROC)  curves  as  a  means  of  identifying  in  a  rigorous  statistical  manner  the  optimal  parameters 
necessary  for  the  classification  scheme  to  operate  efficiently,  including  the  choice  of  alpha  (Type-I 
error)  that  results  in  the  best  trade-off  between  Type-I  and  Type-II  error. 

The  classification  algorithm  can  have  one  of  four  outcomes  for  a  particular  test  case.  A  true 
positive  occurs  when  the  bond  being  interrogated  is  damaged  and  the  algorithm  also  concludes  that 
it  is  damaged.  If  the  algorithm  concludes  that  this  damaged  bond  is  undamaged  that  is  referred  to  as 
a  false  negative.  A  false  positive  occurs  when  the  bond  being  interrogated  is  undamaged  and  the 
algorithm  concludes  that  it  is  damaged.  Similarly,  if  the  algorithm  correctly  concludes  that  this 
undamaged  bond  is  undamaged  that  is  known  as  a  true  negative.  ROC  curves  are  formed  by  plotting 
the  rate  of  true  positives  found  by  the  classification  scheme  (true  positives  divided  by  true  positives 
plus  false  negatives)  against  the  false  positive  rate  (false  positives  divided  by  false  positives  plus 
true  negatives).  Essentially,  the  more  area  under  a  ROC  curve  the  better  job  of  correctly  identifying 
bond  state  that  particular  parameter  is  doing  because  there  will  be  a  high  number  of  true  positives 
compared  to  a  low  number  of  false  positives. 


77 


The  bond  state  assessment  algorithm  is  used  for  an  array  of  signal  parameters  and  the 
overall  ability  of  each  set  of  parameters  to  correctly  identify  each  bond  state  is  then  examined.  In 
this  study  a  number  of  excitation  center  frequencies  from  100  kHz  to  350  kHz  are  considered,  as 
well  as  frequency  ratios  ranging  from  0.15  (slightly  chaotic)  to  0.50  (highly  chaotic).  When 
considering  the  feature  extraction  parameters  a  number  of  AR  model  orders  from  5  to  30  were 


employed.  Figure  40  depicts  the  ability  of  the  VCC/EVS-based  statistical  classification  scheme  to 
correctly  identify  bond  condition  by  plotting  the  receiver  operating  characteristic  (ROC)  curve.  This 
plot  depicts  the  true  positive  rate  against  the  false  positive  rate  for  all  possible  classifications  (both 
undamaged  and  damaged).  Essentially  the  more  area  under  a  ROC  curve  the  better  job  of  correctly 
identifying  bond  state  that  particular  parameter  is  doing  because  there  will  be  a  high  number  of  true 


positives  compared  to  a  low  number  of  false  positives  (Type  I  error). 


False  Positive  Rate 


Figure  40.  ROC  curves  for  (a,  top  left)  Plate  1  center  frequencies,  (b,  top  right)  Plate  2  center 


frequencies,  (c,  lower  left)  Plate  1  frequency  ratio  R  and  (d,  lower  right)  Plate  1  AR  model  order. 
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Figure  40  (top  left)  shows  various  ROC  curves  acquired  by  changing  the  center  frequency  of  the 
excitation  signals  on  Plate  1.  Performance  of  the  bond  state  assessment  algorithm  increases  fi*om  a 
center  frequency  of  100  kHz  up  until  300  kHz,  at  which  point  the  performance  degrades  as  the 
center  frequency  continues  to  increase  to  350  kHz  (only  center  frequencies  from  250  kHz  to  350 
kHz  are  shown  here).  This  means  that  for  the  particular  geometry  of  the  test  specimen  and  with 
MFC  patches  affixed  to  the  same  side  of  the  plate  as  the  bond  line  that  a  center  frequency  of  300 
kHz  is  best  suited  for  bond  state  assessment.  Figure  40  (top  right)  depicts  the  same  change  in  center 
frequencies  as  the  top  left  plot  except  that  the  data  is  from  Plate  2  which  has  the  MFC  patches 
bonded  on  the  opposite  side  of  the  bonded  spar.  There  is  similar  behavior  to  the  previous  plot  in 
which  performance  increases  from  100  kHz  to  a  particular  optimal  frequency  and  then  decreases 
again.  However  in  this  case  the  optimal  center  frequency  is  275  kHz.  This  change  is  due  to  the 
different  way  the  guided  wave  propagates  through  the  bond  when  introduced  to  the  plate  on  the 
opposite  side  of  the  spar. 

Figure  40  (bottom  left)  shows  the  ability  of  the  statistical  classification  scheme  employed 
together  with  chaotic  insonification  to  correctly  identify  bond  condition  for  various  frequency  ratios 
on  Plate  1.  A  frequency  ratio  of  0.3  appears  to  be  best  suited  for  bond  state  assessment.  The  ROC 
curves  for  Plate  2  frequency  ratios  are  very  similar  and  again  show  that  a  frequency  ratio  of  0.3  is 
optimal.  Figure  40  (bottom  right)  depicts  the  effect  of  AR  model  order  on  the  damage  identification 
ability  of  the  statistical  classification  scheme.  All  AR  model  orders  appear  to  work  approximately 
equally  well,  except  for  an  order  of  5,  with  an  AR  model  order  of  20  performing  the  best  among  the 
group.  Again  Plate  2  shows  similar  results  for  the  AR  model  order.  It  is  not  possible  to  create  a 
ROC  curve  for  Plate  3,  as  there  are  no  true  positive  results  because  there  is  no  damage  (and  thus,  no 
hypothesis  test).  Figure  41  shows  the  optimal  ROC  curves  for  Plate  1  and  Plate  2. 
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Figure  41.  Optimal  ROC  curves  for  (a,  left)  Plate  1  and  (b,  right)  Plate  2. 

As  with  the  single  plate,  this  study  again  employs  the  use  of  a  thermal  chamber  to  test  the 
efficacy  of  the  above  outlined  method  for  multiple  plates  for  temperatures  varying  from  -40°C  to 
40°C,  a  range  that  covers  most  flight  and  ground  conditions.  The  optimal  set  of  parameters  that  will 
be  used  in  the  following  experiments  is  a  center  frequency  of  300  kHz  for  Plate  1  and  Plate  3  (MFC 
patches  on  the  same  side  as  the  bond)  and  a  center  frequency  of  275  kHz  for  Plate  2  (MFC  patches 
on  the  opposite  side  as  the  bond).  A  frequency  ratio  of  0.3  will  be  used  for  all  specimens  as  well  as 
an  AR  model  order  of  20.  These  are  the  most  optimal  parameters  that  were  identified  using  ROC 
curves  in  the  previous  section. 

An  alpha  level  of  0.05  (5%  type  I  error)  is  used  to  choose  the  lower  confidence  limit  for  the 
undamaged  VCC  distribution.  Figure  42  shows  the  percentage  of  outliers  for  each  bond  state  being 
inspected  on  Plate  1  and  Plate  2.  This  percentage  of  outliers  is  plotted  against  temperature  and  the 
“training  database”  of  AR  coefficient  vectors  has  been  acquired  at  each  temperature.  If  the 
percentage  of  outliers  is  above  the  5%  type  I  error  that  has  been  specified,  the  bond  state  assessment 
algorithm  has  determined  that  the  unknown  bond  being  inspected  is  damaged. 
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Figure  42.  Percentage  of  outliers  vs.  temperature  for  each  bond  condition  on  (left)  Plate  1,  (right) 
Plate  2. 


Figure  42  shows  that  for  every  temperature  the  undamaged  path  and  paths  with  disbonds  are 
correctly  classified  because  the  blue  undamaged  line  is  always  below  the  decision  boundary  and  the 
red  and  green  disbond  lines  are  always  above  the  decision  boundary.  It  is  therefore  possible  to 
correct  for  temperature  variations  by  storing  baseline  AR  coefficients  for  the  range  of  temperatures 
that  the  intended  apparatus  is  likely  to  see  during  real  world  applications.  However  it  is  not 
necessary  to  store  AR  coefficients  for  all  possible  temperatures  because  only  a  few  temperatures 
need  to  be  stored  and  AR  coefficients  for  intermediate  temperatures  can  be  linearly  interpolated 
between  existing  baseline  AR  coefficients.  This  study  again  uses  the  5tempbase  and  9tempbase 
methods  of  AR  coefficient  interpolation  outlined  previously.  Testing  is  carried  out  at  5  degree 
Celsius  temperature  increments. 

Figure  43  depicts  the  same  results  as  are  seen  Figure  42  except  that  AR  coefficient  values 
are  linearly  interpolated  at  temperatures  for  which  no  baseline  data  is  present.  The  5tempbase 
method  depicted  in  the  left-hand  plot  of  Figure  43  works  well  except  for  two  temperatures  at  which 
the  number  of  outliers  for  the  undamaged  path  rises  above  the  critical  threshold.  The  9tempbase 
method  shown  in  the  right-hand  plot  has  100%  correct  classification  at  all  temperatures.  This  study 
reaffirms  that  the  best  method  for  temperature  correction  is  to  store  baseline  coefficients  for  a  small 
number  of  temperatures  that  span  the  known  operational  values  of  the  structure  and  then  to  linearly 
interpolate,  via  a  lookup  table,  the  AR  coefficients  necessary  to  correctly  classify  bond  condition  for 
any  intermediate  temperatures. 
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Figure  43.  Percentage  of  outliers  for  Plate  2  using  the  Stempbase  (a,  left)  and  9tempbase  (b,  right) 
methods  of  AR  coefficient  interpolation. 


Having  already  determined  what  paths  have  some  form  of  disbond  or  damage,  the  second 
bond  state  assessment  algorithm  employed  in  this  study  seeks  to  classify  the  size  of  the  disbond. 
Figure  44  shows  the  percentage  of  total  correct  votes  for  each  disbond  size  using  either  the 
Stempbase  or  9tempbase  AR  coefficient  interpolation  method  for  Plate  1  and  Plate  2.  The 
Stempbase  and  9tempbase  methods  appear  to  work  equally  well  and  the  large  disbond  is  correctly 
identified  in  almost  100%  of  individual  cases.  The  small  disbond  is  occasionally  misclassified  but 
the  lowest  result  at  1  S®C  on  Plate  2  is  still  greater  than  90%  correct  classification. 
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Figure  44.  Percentage  of  total  correct  votes  for  each  disbond  size  using  the  Stempbase  and 
9tempbase  AR  coefficient  interpolation  methods  on  Plate  1  (a,  left)  and  on  Plate  2  (b,  right). 
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4.3.4.  Other  Sources  of  Variability 


It  is  possible  to  observe  the  resulting  effects  of  within-unit  variability  by  examining  Plate  3, 
in  which  all  paths  were  undamaged.  Figure  45  shows  the  number  of  outliers  using  the  VCC  feature 
for  Plate  3.  In  this  figure,  one  path  was  chosen  to  be  the  baseline  undamaged  path  and  the  other  two 
paths  were  used  as  test  cases.  The  number  of  outliers  for  each  of  the  two  test  undamaged  paths  is 
above  the  5%  threshold  but  are  less  than  7.5%  and  much  less  than  the  number  of  outliers  seen  in  the 
paths  with  disbonds  on  Plate  1  and  Plate  2  shown  Figure  42.  This  result  only  comes  from  one  test 
structure,  but  it  appears  that  it  may  be  possible  to  build  a  safety  factor  (in  this  case  1.5)  into  the 
lower  confidence  limit  to  account  for  the  various  sources  of  within-unit  variability  on  a  particular 
structure  and  still  be  able  identify  when  and  where  the  state  of  the  bond  has  been  changed.  This 
result  is  important  because  it  means  it  may  not  be  necessary  to  store  baseline  data  for  every  MFC 
actuator/sensor  pair  in  a  structure,  as  it  may  be  possible  to  generalize  baseline  data  from  elsewhere 
in  the  structure  as  long  as  the  geometries  are  similar. 


Temperature  (C) 

Figure  45.  Percentage  of  outliers  vs.  temperature  for  each  bond  condition  on  Plate  3. 

Unit-to-unit  variability  arises  due  to  inconsistent  manufacturing  processes  from  specimen  to 
specimen,  which  compound  the  within-unit  sources  of  variability  previously  discussed.  In  well- 
controlled  automated  industrial  applications  with  very  tight  tolerances  or  requirements,  unit-to-unit 
variability  is  typically  lower  than  within-unit  variability.  However,  in  this  study  each  specimen  was 
manufactured  individually  by  hand,  and  consequently  the  unit-to-unit  variability  is  considerably 
larger  than  within-unit  variability.  Figure  46  illustrates  the  large  effect  of  this  manufacturing 
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inconsistency.  The  left-hand  plot  shows  averaged  AR  coefficients  that  are  acquired  using  an  AR 
model  order  of/?=10  for  the  three  different  bond  conditions  on  Plate  1  and  the  undamaged  bond 
condition  on  Plate  3  that  is  in  the  same  geometric  location  as  the  undamaged  bond  condition  on 
Plate  1 .  Tbe  plot  reveals  that  the  AR  coefficient  vector  from  the  undamaged  bond  condition  on  Plate 
3  differs  as  much  from  the  AR  coefficient  vector  from  the  undamaged  bond  condition  on  Plate  1  as 
the  AR  coefficient  vectors  from  the  two  disbond  conditions  on  Plate  1.  This  large  difference 
between  the  AR  coefficient  vectors  from  the  undamaged  bond  condition  in  the  same  geometric 
location  on  each  plate  shows  that  the  unit-to-unit  variability  is  too  large  to  use  baseline  data 
gathered  on  one  specimen  to  discern  the  damage  state  of  a  second  specimen.  TTiis  conclusion  is 
depicted  graphically  in  the  right-hand  plot  of  the  figure,  which  shows  the  percentage  of  outliers  vs. 
temperature  for  each  bond  state  on  Plate  1  using  the  undamaged  bond  condition  on  Plate  3  as  the 
baseline  undamaged  condition.  The  percentage  of  outliers  for  the  undamaged  bond  condition  on 
Plate  1  is  well  above  the  decision  boundary  and  would  therefore  be  classified  as  damaged  if  the 
undamaged  bond  condition  from  Plate  3  were  used  as  a  baseline.  TTius  it  is  concluded  that  for  the 
manufacturing  techniques  used  in  this  study  (individual  test  specimen  construction  by  a  non-expert 
human  technician)  that  unit-to-unit  variability  is  too  large  to  allow  for  baseline  data  gathered  from 
one  plate  to  be  used  to  diagnose  the  bond  condition  of  a  second  plate.  TTiis  result  does  not  preclude 
the  possibility  that  with  more  controlled  manufacturing  processes  that  unit-to-unit  variability  could 
be  small  enough  that  baseline  data  from  only  one  specimen  would  be  needed  to  properly  assess  the 
bond  condition  of  a  specimen  lot. 


p 

u 


AR  Coefficient  Temperature  (C) 

Figure  46.  Averaged  AR  coefficient  vector  comparison  for  each  bond  condition  on  Plate  1  vs. 
undamaged  bond  condition  on  Plate  3  (left).  Percentage  of  outliers  vs.  temperature  for  each  bond 
condition  on  Plate  1  using  baseline  data  on  Plate  3  as  the  undamaged  baseline  (right). 
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4.3.5.  Section  Summary 

This  section  has  shown  that  the  structural  health  monitoring  capability  of  chaotically 
modulated  ultrasonic  waves  that  are  imparted  to  a  composite  bonded  wing  skin-to-spar  structure 
through  a  piezoelectric  patch.  The  classification  damage  detection  scheme  was  shown  to  be 
effective  in  identifying  bolt  preload  configuration  in  experiments  on  single  and  multiple  joint 
structures.  For  the  single  composite  structure  experiment  in  Section  0  the  damage  size/type 
classification  scheme  is  shown  to  be  effective  in  identifying  various  bond  damage  states  including 
different  disbond  sizes  as  well  as  poorly  cured  bonds  and  has  the  ability  to  compensate  for 
temperature  variation.  This  study  has  also  shown  the  ability  of  a  two-part  bond  state  assessment 
algorithm  to  correctly  classify  the  existence  and  size  of  a  disbond  within  a  composite  bonded  wing 
skin-to-spar  structure.  The  use  of  vector  consistency  criterion  (VCC)  values  derived  from 
comparison  of  AR  coefficient  vectors,  in  combination  with  outlier  analysis  based  on  extreme  value 
statistics,  allows  damage  existence  classification  to  be  possible  even  without  baseline  data  from 
every  MFC  sensor/actuator  pair  in  a  structure.  Disbond  size  can  then  be  determined  in  a  supervised 
learning  manner  by  employing  a  classification  technique  derived  from  information  theory.  This 
method  also  has  the  ability  to  compensate  for  temperature  variability  provided  baseline  AR 
coefficients  have  been  recorded  at  a  suitable  number  of  temperatures  that  span  the  operating  range 
of  the  structure  (and  that  temperature  is  independently  measured  in  situ  during  operation,  which  is, 
in  practice,  typically  done). 
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